# Soot Formation Models Assessment in Turbulent Diffusion Jet Flames: A RANS-based Comparison

Sebastian Valencia, Cesar Celis, Luís Fernando Figueira da Silva

► To cite this version:

Sebastian Valencia, Cesar Celis, Luís Fernando Figueira da Silva. Soot Formation Models Assessment in Turbulent Diffusion Jet Flames: A RANS-based Comparison. 2024. hal-04717636

HAL Id: hal-04717636

https://hal.science/hal-04717636v1

Preprint submitted on 2 Oct 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Soot Formation Models Assessment in Turbulent Diffusion Jet Flames:
A RANS-based Comparison

Sebastian Valencia<sup>a*</sup>, Cesar Celis<sup>a</sup>, and Luís Fernando Figueira da Silva<sup>b</sup>

<sup>a</sup>Mechanical Engineering Section, Pontificia Universidad Católica del Perú, Lima, Peru; <sup>b</sup>Institut Pprime, CNRS, ISAE-ENSMA, Université de Poitiers, Futuroscope Chasseneuil, France

*corresponding author. Email: svalenciar@pucp.edu.pe
Soot Formation Models Assessment in Turbulent Diffusion Jet Flames:
A RANS-based Comparison

Due to the intricate interaction between turbulence, chemical kinetics, radiation,
and soot particle dynamics, modelling soot formation processes in flames is a
challenging task. To predict the level of soot formed, it is essential to accurately
capture all stages of soot formation and oxidation. Using a RANS approach, this
study focuses on the implementation, within the computational open-source tool
OpenFOAM, and comparison of three detailed soot formation models, (i) the
Interpolative Closure Method of Moments (MOMIC), (ii) the Hybrid Method of
Moments (HMOM), and (iii) the Discrete Sectional Method (DSM), as well as a
semi-empirical two-equation model. Both the combustion process and the
formation of soot precursors in the gas phase are described using the Steady
Laminar Flamelet model and a detailed chemical kinetic mechanism. Radiation
effects are modelled using the optically thin method. The computational results
obtained here are compared with the experimental data characterizing the
Adelaide ENH1 jet flame and other past numerical results. The results reveal
significant differences in soot formation source terms among the models and for
each of the soot formation stages. DSM best matches the experimental peak soot
position, which is attributed to its modelling of condensation and surface growth
occurring downstream compared to MOMIC and HMOM.

Keywords: RANS; soot modelling; turbulent diffusion flames; method of
moments; discrete sectional method
# 1. Introduction

Over the years, substantial progress has been made in enhancing the comprehension of soot formation. Nevertheless, the intricate process of carbonaceous particle formation within hydrocarbon flames has persisted as a perennially unsolved enigma [1]. Until now numerous aspects of this phenomenon remain shrouded in uncertainty, particularly those related to the transition from gaseous phase species to the inception of the initial soot particles called nuclei. A widely accepted conceptual framework for soot formation and transport involves a process that starts with the conversion of gas-phase hydrocarbons into geometrically compact polycyclic aromatic hydrocarbon (PAH) precursors. These PAH precursors evolve then into solid nuclei that undergo gradual enlargement through the deposition and condensation of PAH, known as surface growth. This gradual growth process leads to the formation of larger fractals referred as aggregates. Ultimately, these relatively large fractals undergo defragmentation and consumption through chemical reactions with oxygen and hydroxyl radicals via oxidation processes [2].

Predicting soot formation in turbulent reacting flows entails numerous complexities stemming from the intricate interplay between turbulence, chemistry, soot formation, and radiation [3]. To characterize soot formation processes, it is imperative to properly model the gas-phase kinetics associated with both soot precursors and oxidizers. Among the widely employed combustion models used for numerically modelling complex turbulent non premixed flames is the Steady Laminar Flamelet Model (SLFM) [4]. In particular, the integration of tabulated chemistry with soot formation models has seen extensive usage in the simulation of sooting flames [5-7]. However, it is important to notice that PAH precursors, responsible for both particle inception and growth, are subject to a slow chemical process, and a longer time is needed to respond to rapid
changes in the reacting flow field. Consequently, the relatively longer timescales of
PAH formation, when compared to other gas-phase species, often violate the
assumption of fast chemistry inherent in flamelet-based models [8]. An alternative
solution to address this limitation is to calculate PAH mass fractions using a transport
equation rather than relying on data from a flamelet database [6,9].

Due to the large number of soot particles, which molar fraction is on the order of ppm, tracking individual particles is a challenge in modeling [10]. Therefore, statistical descriptions of the evolution of soot particle populations are often used, which consist in obtaining statistically the particle size distribution function (PSDF) that describes the concentration of soot particles as a function of their size. Broadly speaking, soot formation models are classified into three subgroups [11]: (i) empirical models, (ii) semi-empirical models, and (iii) detailed models. Detailed soot formation models are the most studied ones in recent years. These models include the Method of Moments (MOM) [9,12,13], the Sectional Method [14–17], and the stochastic Monte Carlo method [18,19]. Among these, the Method of Moments and the Sectional Method are the most often applied in turbulent reacting flows. This occurs because of the relatively high computational cost of the Monte Carlo method, which limits it to canonical 0-D and 1-D configurations, and currently makes it unfeasible for multidimensional turbulent flames. MOM-based soot formation models offer in turn a cost-effective means to describe soot morphology. However, when reconstructing the PSDF, they require a closure model for unsolved higher-order statistical moments, which increases their mathematical complexity [20]. Finally, Discrete Sectional Method (DSM) based soot models enable a discretized representation of the soot PSDF in the particles volume space. Hence, the PSDF can be accurately reproduced. Nevertheless, in order to effectively capture the soot PSDF, DSM-based models demand a relatively large
number of sections. This can increase substantially the computational requirements, particularly when dealing with multi-dimensional cases [20].

Accordingly, the main goal of this work is to assess, using a RANS/flamelet approach, the capabilities of four different soot formation models, coupled to an optically thin radiation model (as soot is as a strong emitter and absorber of radiative heat [2]), to predict the soot formed in turbulent diffusion flames. The first soot formation model used here is the classical semi-empirical two-equation (2EQ) one [21], which is based on acetylene as soot precursor. The second model in turn is the Method of Moments with Interpolative Closure (MOMIC) [12], which, to obtain higher order statistical moments and reconstruct the PSDF, uses a population balance approach that solves a set of moment equations with an interpolative closure scheme. The third model used here, the Hybrid Method of Moments (HMOM) [9], which is based on both the MOMIC model and the Discrete Quadrature Method of Moments (DQMOM) [22] one, combines the accuracy of the former with the robustness of the latter. This is a bivariate soot model based on soot particles surface and volume related parameters. The fourth and last model studied in this work is the Discrete Sectional Method (DSM) [16], which divides the particle size distribution into a finite number of discrete sections and solves transport equations for each of these sections. It is worth noticing that, regardless of the soot model employed, the soot source terms exhibit a close dependency on the soot precursor employed within each model. Due to its relatively lower computational cost when compared to more computationally expensive techniques such as LES, a RANS approach has been preferred in this work. This is supported by the fact that, in the present context of model development or implementation, which requires tens or even hundreds of numerical simulations, the computational cost is a critical issue.
In this work, the referred four models are applied to model a turbulent diffusion sooting
flame (Adelaide ENH1 jet flame) [23] subject of study of the International Sooting
Flame (ISF) workshop. The ISF workshop describes indeed several efforts carried out in
the past to model this flame, where each modelling approach has its own strengths and
weaknesses. Notice that the validation process of the numerical results obtained in this
work encompasses an examination of the associated fluid flow dynamics through the
scrutiny of velocity and velocity fluctuation profiles within the jet discharge zone, as
well as the axial temperature profile along the centerline. Furthermore, a comparative
analysis of the results obtained here is also conducted using other numerical results
obtained in the past employing a LES/Radiation Flamelet Progress Variable
(RFPV)/HMOM approach [24]. These last comparisons highlight the capabilities of the
(RANS) model employed in this work to accurately characterize the studied jet flame.
Regarding soot predictions, the ability with which these soot formation models predict
the peak soot concentration is specifically addressed.

The remaining of this work is organized as follows. In Section 2, the main soot
formation and reduction mechanisms are discussed. The mathematical models
employed here, including those related to turbulence, combustion, and radiation, are
detailed in Section 3. In Section 4, the four soot formation models studied in this work
are described in detail. The numerical model, including the flame configuration
accounted for, is highlighted in Section 5. In Section 6, the numerical results obtained
here are compared against experimental data and LES results. Finally, in Section 7, the
main conclusions of the preset work are summarized.

## 2. Soot formation and reduction mechanisms

There is a widespread consensus on the existence of six physical/chemical processes
directly involved in the formation and oxidation of soot, namely (1) pyrolysis, (2) inception, (3) coalescence, (4) surface growth, (5) aggregation, and (6) oxidation. Due to the occurrence of chemical reactions and the increase in the reacting flow temperature, the first stage, called pyrolysis, involves the decomposition of fuel into small chains of HC radicals [25]. These HC radicals promote the formation of acetylene ($C_2H_2$) and polycyclic aromatic hydrocarbons (PAH) which, according to several authors, are possibly the main responsible for the formation of soot [12]. The formation of PAH begins when the first aromatic ring, benzene ($C_6H_6$ or A1) appears. Over the years, different pathways for the formation and growth of PAH have been proposed. One of the most universally accepted pathways is the *H*-Abstraction-$C_2H_2$-Adition (HACA) mechanism [26]. The HACA mechanism consists of two chemical reactions and involves the abstraction of a hydrogen (*H*) atom and the addition of acetylene ($C_2H_2$). The first reaction is the abstraction of an *H* radical from a precondensed aromatic ring and the change of the aromatic ring to a more active polycyclic hydrocarbon radical. In the second reaction, an acetylene molecule is added to the site of the free radical of the ring or aromatic compound.

The coalescence of PAHs is the most prevalent inception mechanism for soot formation. According to this mechanism, through random Brownian collision, large PAHs collide and coalesce to create the first and smallest soot particles [26]. The elements formed due to the coalescence of PAHs are called dimers. Dimers are solid-phase elements formed by stacking two PAH monomers. Through experimental studies [27], it has been shown that PAH dimers exist in the proximity of the inception zone of a laminar methane diffusion flame. In a more recent study [28], it was concluded that soot inception is also produced through a mechanism like HACA. Unlike other inception mechanisms, this study showed that a more recurrent activation growth favors the
efficiency of PAH growth processes. However, to the best of the authors knowledge there is still no computational model that accounts for this new inception mechanism. In this work, an inception model involving the collision of PAHs is used for MOMIC and DSM, whereas a PAH dimerization model [29] is used in HMOM. Once soot is formed, it evolves due to physical interaction and chemical reactions with species in the gas phase. Physical processes include PAH coagulation and condensation, whereas chemical processes include surface reactions and oxidation.

The condensation mechanism involves the increase in size of soot particles due to their collision with PAHs or dimers. When a PAH or dimer collides with a soot particle, it is most likely that the weakly bonded dimer breaks apart, and the PAH molecules are deposited on the soot particle surface [30]. As a result, the volume of the soot particle participating in the collision increases, and it becomes more spherical due to the mass condensing on its surface. In this work, the condensation phenomenon is accounted for in HMOM and DSM only.

Surface growth describes chemical reactions occurring on the surface of a soot particle due to the presence of other chemical species available in the reacting environment. This surface growth is primarily driven by the HACA mechanism [26], in which the C-H bonds on the soot particle surface provide available reaction sites for absorbing gaseous species, particularly acetylene. Acetylene ($C_2H_2$) is considered indeed as one of the most important gaseous species in this soot formation mechanism [31].

The coagulation process involves in turn the collision of soot particles. When the resulting soot particle from the collision of two smaller particles is spherical, this coagulation is typically referred to as coalescence. However, when the resulting soot particle preserves the essential structure of the colliding particles and forms an aggregate, the coagulation is known as aggregation. The coagulation rate depends on
the size and shape of the colliding particles [9]. Thus, to model the coagulation regime,
it is necessary to estimate the soot particles Knudsen number. In this work, the
aggregation phenomenon was exclusively incorporated into HMOM, as it is a bivariate
model.

Finally, soot oxidation is a process that reduces the mass of soot particles, through the
partial or complete conversion of solid soot particles into gaseous-phase products such
as CO and CO2. This reduction in the size of soot particles occurs through surface
oxidation originated from the collision of these particles with OH radicals and O2
molecules [32]. It is worth noticing here that the oxidation of soot particles can lead to
their fragmentation. The fragmentation process further decreases the size of soot
particles as they break up into multiple smaller particles. Nevertheless, due to the
complexity of the processes involved in soot formation, most current numerical models
do not account for soot particles fragmentation phenomena [33].

**3. Mathematical Modeling**

The main details of the mathematical model employed in this work to describe turbulent
reacting flows featuring soot formation are presented in this section. Specifically, the
different modeling approaches used to model turbulence, chemical reactions, and
radiation are briefly described.

**3.1. Turbulence Modeling**

This work focuses on the analysis of low Mach number turbulent reacting flows
characterized by variable density. To this end, as turbulence modeling approach, the
Reynolds Averaged Navier-Stokes (RANS) one is utilized here, in conjunction with the
standard *k*–ε turbulence model [34]. This turbulence modeling approach is based on the
temporal averaging of the flow governing equations for mass, momentum, and energy, with two additional transport equations solved to evaluate the turbulent kinetic energy (k) and the turbulent dissipation rate (ε). In order to improve the prediction of the axial velocity decay rate of a round jet, it has been suggested [35] that the coefficient C₁ (present in the transport equation for ε [34]) should be modified, for instance by selecting a value of 1.6 [36]. For the purposes of this work, a value of 1.53 has been chosen, as larger values were found to result in soot leaving the flame tip.

## 3.2. Combustion Modeling

The combustion model employed here is the steady laminar flamelet model SLFM [37]. This particular model accounts for a turbulent flame as an ensemble of thin, one-dimensional flamelets within the turbulent flow field. The SLFM relies on a library of flamelets, which contains pre-computed solutions of the steady, one-dimensional flame equations at various mixture fraction and scalar dissipation rate values. The mixture fraction, Z, represents indeed the local fuel-to-oxidizer ratio in the combustion process. Thus, the thermodynamic properties and mass fractions of each chemical species are extracted from the flamelet library, which is constructed as a function of the mixture fraction, Z, and the scalar dissipation rate, χ. These two parameters are obtained by solving the steady flamelet equations [4],

$$ -\rho \frac{\chi}{2} \frac{\partial^2 \phi}{\partial Z^2} = \omega, \qquad (1) $$

where φ is the vector containing the mass fractions of the chemical species and the temperature, whereas ω denotes their corresponding source terms. The scalar dissipation rate term χ appears as an external parameter and is defined as [4],

$$ \chi = 2D_z(\nabla Z)^2, \qquad (2) $$
where $D_z$ is the mixture fraction molecular diffusivity. The scalar dissipation rate $\chi$ is parametrized as a function of its value at stoichiometry, $\chi_{st}$. Therefore, in the flamelet library, the temperature and mass fractions of the species are tabulated as a function of $Z$ and $\chi_{st}$. The average values of the vector $\phi$ are computed using a joint probability density function (PDF), $P(Z, \chi_{st})$. In this work, statistical independence is assumed, $P(Z, \chi_{st}) = P(Z)P(\chi_{st})$. For the scalar dissipation rate, the form of the PDF is assumed to be a Dirac delta function, whereas for the mixture fraction a beta-type one.

Finally, notice that, to use the SLFM, it is necessary to solve two additional transport equations, (i) one for the mean mixture fraction $\tilde{Z}$, and (ii) another for its variance $\tilde{Z}''^2$. Since $\chi$ is obtained from $Z$ (Eq. (2)), there is no need to solve an additional transport equation for this parameter. The stoichiometric scalar dissipation rate is in turn obtained using a classical similarity hypothesis [4],

$$ \chi_{st} = C_{\chi} \frac{\epsilon}{k} \tilde{Z}''^2. \qquad (3) $$

### 3.3. Radiation Modeling

The coupling between the flamelet library and the radiative heat transfer is carried out here using an additional parameter called the enthalpy defect, $H_{def}$ [38],

$$ H_{def} = H - [H_{ox} + Z(H_{fuel} - H_{ox})], \qquad (4) $$

where $H$ stands for the actual enthalpy of the turbulent reacting flow, whereas $H_{ox}$ and $H_{fuel}$ represent the enthalpy of the oxidizer and fuel, respectively. Thus, a transport equation for enthalpy is solved, and the average values of the vector $\phi$ are obtained through a PDF as a function of $Z$, $\chi_{st}$ and $H_{def}$. In this work, the assumed form of the marginal PDF for the enthalpy defect $P(H_{def} | Z, \chi_{st})$ is also a Dirac delta function.
Finally, the source term present in the enthalpy transport equation is computed using the optically thin approximation [39]. This simplified approach accounts for radiation emission only and neglects reabsorption. This method allows carrying out reasonably accurate predictions in various laboratory flames, but overestimates slightly radiation [40]. Specifically, the effect of radiation is introduced through a source term of the form [40],

$$Q_r = -4\sigma a(T^4 - T_{env}^4), \qquad (5)$$

where $\sigma$ is the Stefan-Boltzmann constant, $a$ is the Planck mean absorption coefficient of the gas mixture, and $T_{env}$ is the ambient temperature. In this work, the individual contributions to the radiative heat transfer of $H_2O$, $CO_2$, $CO$, $CH_4$ and soot are accounted for [41].

## 4. Soot modeling in turbulent combustion

Four different soot formation models are studied in this work, (i) the classical semi- empirical two-equation model (2EQ), (ii) the Method of Moments with Interpolative Closure (MOMIC), (iii) the Hybrid Method of Moments (HMOM), and (iv) the Discrete Sectional Method (DSM). Only the main features of each of these four models are highlighted in this section. For a complete model description then, the interested reader should refer to the cited references.

### 4.1. Semi-Empirical Two Equation Model (2EQ)

The two-equation model accounted for here is a semi-empirical one that, to carry out soot predictions, solves two transport equations. It assumes indeed that the soot particles are uniformly distributed and that the formation and oxidation of soot are controlled by the gas-phase reaction rates. More specifically, soot is computed from the solution of
the transport equations for the soot mass fraction $Y_s$, Eq. (6), and the normalized concentration of soot nuclei $b_{nuc}$, Eq. (7),

$$ \frac{\partial \rho Y_s}{\partial t} + \frac{\partial \rho u_j Y_s}{\partial x_j} = \frac{\partial}{\partial x_j} \left( 0.55 \frac{\mu_{eff}}{T} \frac{\partial T}{\partial x_j} Y_s \right) + \dot{M}, \quad (6) $$

$$ \frac{\partial \rho b_{nuc}}{\partial t} + \frac{\partial \rho u_j b_{nuc}}{\partial x_j} = \frac{\partial}{\partial x_j} \left( 0.55 \frac{\mu_{eff}}{T} \frac{\partial T}{\partial x_j} b_{nuc} \right) + \dot{N}. \quad (7) $$

Since soot is not in the gas phase, molecular diffusion is not included in this model. Instead, thermophoresis, the first term on the right-hand side of Eqs. (6) and (7), is accounted for. In these equations $\mu_{eff}$ stands for the effective viscosity and the concentration rate of nuclei $\dot{N}$ and the mass fraction rate $\dot{M}$ are source terms used to model the contribution to soot formation of soot inception, coagulation, surface growth, and oxidation. These source terms are detailed in reference [21]. The 2EQ formulation studied here incorporates model constants aimed at achieving a satisfactory fit for the soot volume fraction distribution. It is worth noticing here that no modifications were made to these constants in this work [21].

## 4.2. Aerosol Dynamics

MOM and DSM based models rely on the transport of the soot particle size distribution function (PSDF). Following this approach, the temporal evolution of the ensemble of soot particles is described by the population balance equation (PBE), originally based on the discrete coagulation equations [42],

$$ \frac{\partial N_1}{\partial t} = - \sum_{j=1}^{i-1} \beta_{1,j} N_1 N_j, \quad (8) $$

$$ \frac{\partial N_i}{\partial t} = \frac{1}{2} \sum_{j=1}^{i-1} \beta_{j,i-j} N_i N_{i-j} - \sum_{j=1}^{\infty} \beta_{i,j} N_i N_j. \quad (9) $$
In these two equations, $N_i$ denotes the number density of soot particles of size $i$, and $\beta_{i,j}$ is the collision frequency between particles of size $i$ and $j$. The right-hand side positive terms in Eqs. (8) and (9) correspond to the particles added to the class of size $i$, whereas the negative terms refer to the loss of particles from this class. Due to particle collision, these last particles are transferred to a larger size class. The collision coefficient $\beta_{i,j}$ depends, in turn, on the coagulation regime, classified according to the Knudsen number $k_n$ [12]. Depending on the value of $k_n$, particle collision can be in three different regimes. Indeed, when the Knudsen number is greater than 10, the regime is free molecular, whereas for Knudsen numbers less than 0.1, the regime is continuous. For intermediate Knudsen numbers, a transient regime, in which the collision frequency factor is approximated as the harmonic mean of the other two regimes, is accounted for.

### 4.3. Methods of Moments with Interpolative Closure (MOMIC)

In this soot formation model, the PSDF is characterized exclusively by a finite number of low order statistical moments. It should be noticed however that this model is subject to some limitations. For instance, it is unable to account for the bimodality of the PSDF [30]. In this approach, the rth-order moment of the number density of particles $N_i$ along with the corresponding reduced moment $\mu_r$ are described according [43],

$$M_r = \sum_{i=1}^{\infty} m_i^r N_i, \qquad (10)$$

$$\mu_r = \frac{M_r}{M_0}, \qquad (11)$$

where $m_i$ and $N_i$ represent the mass and number density of soot particles of size $i$, respectively. Both the zeroth-order moment $M_0$ and the first-order moment $M_1$ have physical interpretations in this case. Indeed, $M_0$ is the total number density of particles,
whereas $M_1$, when divided by the soot density, represents the soot volume fraction. In this model, the transport equation for the rth-order moment is given by [12],

$$ \frac{\partial \rho M_r}{\partial t} + \frac{\partial \rho u_j M_r}{\partial x_j} = \frac{\partial}{\partial x_j} \left( 0.55 \frac{\mu_{eff}}{T} \frac{\partial T}{\partial x_j} M_r \right) + \dot{M}_r, \quad (12) $$

where thermophoresis is accounted for, and $\dot{M}_r$ represents the source term and includes the contribution from soot particle inception, coagulation, surface growth and oxidation.

A general form for these source terms is expressed as [12],

$$ \dot{M}_r = \dot{M}_r^{nuc} + \dot{M}_r^{coag} + \dot{M}_r^{sg} + \dot{M}_r^{ox}, \quad (13) $$

where $\dot{M}_r^{nuc}$, $\dot{M}_r^{coag}$, $\dot{M}_r^{sg}$, and $\dot{M}_r^{ox}$ are the source terms due to inception, coagulation, surface growth, and oxidation, respectively.

The formulation for inception utilized in this work is [26],

$$ \dot{M}_0^{nuc} = \gamma \epsilon \sqrt{\frac{4\pi k_B}{m_C N C_{PAH}}} (d_{PAH} N_{PAH})^2, \quad (14) $$

$$ \dot{M}_r^{nuc} = 2NC_{PAH}\dot{M}_{r-1}^{nuc}, \quad (15) $$

where $m_C$ represents the mass of a carbon atom, $NC_{PAH}$ is the number of carbon atoms in a PAH molecule, $d_{PAH}$ is the diameter of the PAH molecule, $N_{PAH}$ is the number density of PAH molecules, $k_B$ is the Boltzmann constant, $\epsilon$ is the Van der Waals enhancement factor, and $\gamma$ is a sticking coefficient, introduced to reduce both the efficiency of the collisions between small PAH molecules and the soot concentration [29]. As detailed in Section 6, the value of $\gamma$ has been adjusted here to obtain the same maximum peak of soot volume fraction as the experimental data.

The coagulation source terms are in turn described by the expressions given in [12], which depend on the Knudsen number $k_n$. To determine the coagulation source terms indeed, it is necessary to know reduced statistical moments of fractional and negative
order. These reduced statistical moments are determined here from the transported low-order first moments using a Lagrange interpolation method. Good interpolation accuracy has been obtained using three statistical moments only [44]. Therefore, in this work, the number of statistical moments transported is limited to three, M0, M1 and M2, and the missing reduced moments are calculated from,

$$log(\mu_p) = L_p(log(\mu_0), log(\mu_1), log(\mu_2)), \quad (16)$$

where $L_p$ is the Lagrange interpolation function obtained from the first three statistical moments transported here, and $p$ represents the order of the reduced moment, $\mu_p$, to be determined. Surface growth and oxidation are modeled using the HACA mechanism [26]. In this mechanism, reaction rates for both surface growth and oxidation are calculated using Arrhenius equations [45]. Specifically, the rate of soot surface growth due to $C_2H_2$ and the rate of oxidation due to $O_2$ and $OH$ are obtained using the expressions provided in reference [12].

### 4.4. Hybrid Method of Moments (HMOM)

To address the limitation of the MOMIC model, that considers spherical particles only, a modification has been introduced into the associated model formulation [30], such that the moments are represented as functions of both the soot surface and volume, as opposed to the original model formulation that depends on soot mass only. Furthermore, to capture the bimodality of the PSDF, the Direct Quadrature Method of Moments (DQMOM) has been developed [22], in which the PSDF is approximated using a series of Dirac delta functions. However, this last method requires the inversion of a linear system of equations that in many cases is ill-conditioned.

The Hybrid Method of Moments (HMOM) [9] combines thus the accuracy of DQMOM with the simplicity of MOMIC. HMOM uses indeed the DQMOM approach for small
spherical particles and the MOMIC one for large fractal aggregates. To include aggregate chains, the HMOM model describes soot in terms of its volume $V$ and surface area $S$, where the diameter of a primary particle $d_p$ and the number of primary particles in an aggregate $n_p$ are given by [9],

$$d_p = 6VS^{-1}, \qquad (17)$$

$$n_p = \frac{1}{36\pi} V^{-2} S^{-1}. \qquad (18)$$

With this formulation, the statistical moments of the PSDF have two indices, $x$ for the order of the moment in volume, and $y$ for the order of the moment in surface area [9],

$$M_{x,y} = \sum_{i=0}^{\infty} V_i^x S_i^y N_i. \qquad (19)$$

Thus, the transport equation for an $x, y$-th statistical moment $M_{x,y}$ is given by [9],

$$\frac{\partial M_{x,y}}{\partial t} + \frac{\partial u_j M_{x,y}}{\partial x_j} = \frac{\partial}{\partial x_j} \left( 0.55 \frac{\mu_{eff}}{T} \frac{\partial T}{\partial x_j} M_{x,y} \right) + \dot{M}_{x,y}, \qquad (20)$$

where $\dot{M}_{x,y}$ is the source term including inception, condensation, coagulation, surface growth and oxidation. As in the MOMIC model, $\dot{M}_{x,y}$ depends on many moments of different orders, so a Lagrangian interpolation is necessary. Using only three statistical moments, a first-order interpolation polynomial is employed in this work,

$$M_{x,y}^{HMOM} = M_{x,y}^{DQMOM} + M_{x,y}^{MOMIC} = N_0 V_0^x S_0^y + N_L V_L^x S_L^y, \qquad (21)$$

where $N_0$ is the number density of the first nucleated particles, with its respective surface $S_0$ and volume $V_0$ coordinates, highlighting the fact that the size of the first particle is fixed. Additionally, $N_L$ is the number density of large particles, $V_L$ is the average volume of the large particles, and $S_L$ is the average surface area of the large particles. Notice that $N_0$ requires an additional transport equation analogous to Eq. (20).
In the inception model used in HMOM [29], the first soot particles are formed from the collision of two PAH dimers. Thus, the inception source term is expressed as follows [9],

$$ \dot{M}_{x,y}^{nucl} = \frac{1}{2} \beta_N [DIMER]^2 V_0^x S_0^y, \qquad (22) $$

where [DIMER] and $\beta_N$ is represent, respectively, the molar concentration and collision rate of PAH dimers [29]. To obtain the value of [DIMER], it is necessary to determine the dimerization rate $\dot{\omega}_{DIMER}$, which is defined by [29],

$$ \dot{\omega}_{DIMER} = \gamma \sqrt{\frac{4\pi kT}{m_{PAH}}} \left(\frac{6m_{PAH}}{\pi\rho_s}\right)^2 [PAH]^2, \qquad (23) $$

where $m_{PAH}$ is the mass of the PAH and [PAH] its molar concentration. Like in MOMIC, a sticking coefficient $\gamma$ is also introduced here to reduce the efficiency of the collisions between small PAH, to reduce the concentration of DIMERS, and to adjust the peak value of soot volume fraction.

Condensation in turn occurs as a result of the collision between a soot particle and a PAH dimer. Accordingly, the condensation source term is given by [9],

$$ \dot{M}_{x,y}^{cond} = \sum_{i=0}^{\infty} \beta_{C_i} \left( x \frac{\delta V}{V_i} + y \frac{\delta S}{S_i} \right) V_i^x S_i^y [DIMER] N_i, \qquad (24) $$

where $\delta V$ and $\delta S$ are, respectively, the changes in volume and surface area due to condensation. $\beta_{C_i}$ is the collision rate between dimers and soot particles. Given that [DIMER] is already influenced by the sticking coefficient $\gamma$, this coefficient also affects the source term due to condensation.

In addition, the source terms for surface growth are obtained from [9],

$$ \dot{M}_{x,y}^{sg} = \sum_{i=0}^{\infty} k_{sg} \chi S_i ((V_i + \delta V)^x + (S_i + \delta S)^y - V_j^x S_j^y) N_i, \qquad (25) $$
where $\chi$ is the number of active sites per unit area, and $k_{sg}$ is the reaction rate related to the soot surface growth. And, finally, the source term for oxidation is determined by [9],

$$ \dot{M}_{x,y}^{ox} = -k_{ox}\chi \frac{\delta V}{V_0} V_0^x S_0^{y+1} N_0 - \sum_{i=0}^{\infty} k_{ox}\chi \left( x \frac{\delta V}{V_i} + y \frac{\delta S}{S_i} \right) V_i^x S_i^{y+1} N_i, \quad (26) $$

where $k_{ox}$ is the reaction rate related to oxidation. The coagulation source term is defined elsewhere [9].

### 4.5. Discrete sectional method (DSM)

In the discrete sectional soot formation model [16], the PSDF is discretized into a finite number of sections, and, to obtain the final PSDF, a transport equation is solved for each of these sections. In this work, 20 sections were accounted for, where each section $i$ represents a set of particles whose volume is between the volume of the smallest soot particle $v_{min}$ and the volume of the largest particle $v_{max}$. The minimum volume $v_{min}$ is defined here as being equal to the carbon equivalent volume of the number of carbon atoms of the soot precursors (PAH). This carbon equivalent volume is defined for two carbon atoms as,

$$ v_{C_2} = \frac{2MW_C}{N_A\rho_s}. \quad (27) $$

To divide the domain expressed in terms of soot volume, a logarithmic distribution is used in this work [46],

$$ v_i = v_{min} \left( \frac{v_{max}}{v_{min}} \right)^{\frac{i}{n_{bins}}}. \quad (28) $$

Since the sections are continuous, the lower limit of section $i$ is the same as the upper limit of section $i-1$, so section $i$ lies between $v_{i-1}$ and $v_i$. Furthermore, once the discretization of the sections varies linearly in logarithmic scale, a mean particle
diameter takes the form [14],

$$ \bar{v}_i = \frac{v_i - v_{i-1}}{\ln(v_i) - \ln(v_{i-1})}. \qquad (29) $$

Notice that, in each section, the property transported is the soot mass fraction $Y_i$ [14],

$$ \frac{\partial \rho Y_i}{\partial t} + \frac{\partial \rho u_j Y_i}{\partial x_j} = \frac{\partial}{\partial x_j} \left( 0.55 \frac{\mu_{eff}}{T} \frac{\partial T}{\partial x_j} Y_i \right) + \rho_s \dot{Q}_i, \qquad (30) $$

where $\dot{Q}_i$ is the rate of soot production/destruction due to inception, condensation, coagulation, surface growth, and oxidation.

Inception only occurs in the first section and its source term is modelled according to,

$$ \dot{Q}_{1,nuc} = 2\gamma v_{PAH} \beta_{PAH} N_{PAH}^2, \qquad (31) $$

where $\beta_{PAH}$ is the collision rate of two PAHs [16], $N_{PAH}$ is the total number of PAH considered in inception, and $v_{PAH}$ is the total volume of the two colliding PAH. As in the other two detailed soot models studied in this work, the collision of two PAHs depends on the sticking efficiency $\gamma$. The value of $\gamma$ is adjusted in this study to numerically reproduce the maximum soot volume fraction obtained in the experiments. The condensation source term is in turn defined as [14],

$$ \dot{Q}_{i,cond} = \gamma_{cond} v_{PAH} \beta_{PAH} N_{PAH} \int_{v_i}^{v_{i+1}} \beta_{i,PAH} n_i(v) dv, \qquad (32) $$

where $\gamma_{cond}$ is the condensation collision efficiency, $\beta_{i,PAH}$ is the collision frequency function, and $n_i(v)$ is the presumed profile of the soot number density. The coagulation rate is described here in the pure coalescence limit, i.e., it does not account for the aggregation phenomenon. The corresponding source term is described in reference [46]. Surface reactions also occur according to the HACA mechanism [26]. The associated reaction rates are thus computed using equations of Arrhenius form, where the coefficients rates are obtained from [45],
$$ \Delta Q_{i,sg} = 2 v_c k_4 [C_2H_2] \alpha \chi_{C_{soot}*} S_i, \qquad (33) $$

$$ \Delta Q_{i,ox} = 2 m_c k_5 [O_2] \alpha \chi_{C_{soot}*} S_i + m_c \gamma_{OH} [OH] \beta_{OH} S_i N_A, \qquad (34) $$

where $k_4$ and $k_5$ are reaction rates related to surface growth and oxidation, respectively. In addition, $S_i$ is the surface density of soot particles, $\alpha$ is the steric factor, $\chi_{C_{soot}*}$ is the number of radical sites per unit surface area, $\gamma_{OH}$ is the collision efficiency of $OH$, and $\beta_{OH}$ is the collision frequency of $OH$ and soot.

## 4.6. Soot modeling tuned parameters

Detailed soot models aim to eliminate the need for tuned parameters by employing rigorous physics-based methods that maintain accuracy across a broad spectrum of conditions [47]. However, tuning coefficients are still present in all models, especially for soot inception and soot surface mechanisms (condensation, surface growth, and oxidation). For instance, the sticking coefficient $\gamma$, used in inception modeling, is essential in many studies to accurately predict experimental data, so sometimes this parameter is obtained by matching experimental data [48]. Accordingly, the HMOM and MOMIC models studied in this work employs a sticking coefficient obtained from Blanquart and Pitsch [29], where its value depends on the PAH precursor, i.e., 0.002 for naphthalene to 0.025 for pyrene. In this study, since benzene featuring a much higher concentration than naphthalene and pyrene is accounted for as the soot precursor, a much lower sticking coefficient is expected.

Additionally, condensation, surface growth, and oxidation models also include usually tuned factors that influence their final form. For instance, in surface growth models, it is common to adjust the reaction rate $k_4$ (denoted as $k_{sg}$ in HMOM) [45,47]. Similarly, the steric factor $\alpha$ in oxidation models often undergoes modifications. Accordingly,
tuned parameters have also been employed in this work. Table 1 summarizes thus the values of the tuned parameters utilized here. A discussion about the implications of the usage of the referred tuned parameters is provided in Section 6.5.

<table><caption>Table 1. Tuned parameters accounted for here.</caption><thead><tr><th>Model</th><th>Inception</th><th>Condensation</th><th>Surface Growth</th><th>Oxidation</th></tr></thead><tbody><tr><td>MOMIC</td><td>γ = 25×10<sup>-3</sup></td><td>No model</td><td>No modification</td><td>No modification</td></tr><tr><td>HMOM</td><td>γ = 160×10<sup>-6</sup></td><td>Included in [DIMER]</td><td>k<sub>sg</sub> = k<sub>sg</sub> × 0.002</td><td>k<sub>ox</sub> = k<sub>ox</sub> × 0.1</td></tr><tr><td>DSM</td><td>γ = 40×10<sup>-6</sup></td><td>γ<sub>cond</sub> = 10 × 10<sup>-3</sup></td><td>k<sub>4</sub> = k<sub>4</sub> × 0.014</td><td>k<sub>5</sub> = k<sub>5</sub> × 0.02</td></tr></tbody></table>

# 5. Numerical Modelling

The numerical model employed in this work is briefly described in this section. A particular emphasis is put on the solver and the numerical schemes utilized, as well as on the flow configuration and the computational mesh used to obtain the numerical results discussed in this work.

## 5.1. Solver and numerical schemes

The combustion solver used here is an OpenFOAM flamelet solver [49], coupled to the soot formation models developed in this work. This combustion solver is based on the SLFM theory, which, in a context of steady-state RANS simulations, is implemented in OpenFOAM using the SIMPLE algorithm [50]. In this solver, a meticulous integration of all soot models was conducted, entailing the inclusion of the transport equations tailored to each specific term present in the soot models accounted for. Additionally, the model coupling with radiation was established by accounting for the radiative effects of soot, which were calculated within each respective soot model employing Eq. (6).
The library of flamelets required by the solver is constructed in a preprocessing step, using four independent variables, (i) the mixture fraction, (ii) its variance, (iii) the scalar dissipation rate, and (iv) the enthalpy defect. These flamelets are generated by solving the equations describing steady-state laminar counterflow diffusion flames. The generated flamelets cover a wide range of scalar dissipation rates, from equilibrium to extinction. To account for the effects of turbulence on the flamelets laminar profiles, a presumed probability density function (PDF) type approach is employed. Thus, in this combustion solver, the mean mixture fraction, the variance of the mixture fraction, and the mean enthalpy are solved by employing transport equations. In turn, the scalar dissipation field is estimated from both the turbulence field and the mixture fraction variance. Regarding numerical schemes, spatial derivatives are computed using second order Gaussian integration with a flux limiter to avoid unphysical oscillations, and the interpolation scheme employed was linear.

## 5.2. Case study and Computational mesh

The turbulent jet flame Adelaide ENH1 analyzed in this work is subject of study of the ISF workshop [23]. This flame, which is depicted in Figure 1a, utilizes a fuel mixture comprising ethylene, hydrogen, and nitrogen (EHN) with a Reynolds number of 15,000 at the jet exit. The molar concentration of the fuel mixture consists of 39.1% ethylene, 41% hydrogen, and 19.9% nitrogen. The density of the fuel mixture is 0.717 kg/m³, and its dynamic viscosity is 1.215 x 10⁻⁵ kg/m.s. The average velocity of the jet is 56.8m/s, surrounded by an air stream co-flowing with a mean velocity of 1.1 m/s and a turbulent intensity of 1.5%. The jet tube has a conical shape as illustrated in Figure 1b, that has been accounted for in the modelling. Both the fuel mixture and the co-flowing air are at 294 K and 1 bar. The geometrical configuration under study is a straight tube
burner manufactured using aluminum, with an inner diameter of 4.4 mm, wall thickness of 1 mm, and a length of 385 mm. The co-flow duct is a square section of 150x150 mm, and the jet tube exit is raised 18 mm above the square duct.

<figure><img src="image_2.png" /><figcaption>Figure 1 Studied flame and burner: a) Photograph of the jet flame ENH1 [23]. b) Tapered section of pipe jet burner [51].</figcaption></figure>

Experimental measurements available for this flame are velocity and velocity fluctuations located at 2 mm above the burner, axial temperature profile, and profiles of soot volume fraction at different flame locations. Temperature measurements were carried out using a R-type 1 mm diameter thermocouple with a maximum uncertainty of 4%, while laser-induced incandescence (LII) was used to measure the soot volume fraction, with uncertainties of about 25%.

The computational mesh used in this work was generated using the blockMesh utility within OpenFOAM. As illustrated in Figure 2, for the jet-type burner configuration, a wedge-shaped mesh with an internal angle of 5° was constructed. Wedge-type boundary conditions were applied to both faces of the wedge. Fixed values were set at both fuel
and air inlets to maintain temperatures at 294 K and pressures at 101.325 kPa. The fuel inlet velocity is 56.8 m/s, whereas the co-flow air inlet is 1.1 m/s. Turbulent boundary conditions ($k-\epsilon$) were also imposed accounting for turbulent intensities of 5% for fuel and 1.5% for air. More specifically, at the fuel inlet, the turbulent kinetic energy imposed is $11.149 \ 1/m^2s^{-2}$ and the turbulent dissipation rate is $19679.17 \ 1/m^2s^{-3}$. At the air inlet in turn, the turbulent kinetic energy is $0.00041 \ 1/m^2s^{-2}$ and the turbulent dissipation rate is $0.000071 \ 1/m^2s^{-3}$. At the outlet, a zero gradient for exiting flow was prescribed and a patch-normal flux for entering flow was adopted, furthermore, pressure was set to total static pressure, and the remaining fields gradients were set to zero. At the burner wall a no-slip condition was adopted, with wall models used for turbulent kinetic energy and turbulent kinetic dissipation rate. Mean $\tilde{Z}$ was set at 1 for the fuel inlet and 0 for the air inlet, with zero $\tilde{Z}''$ at both inlets. In addition, zero-gradient conditions were applied to all walls for $\tilde{Z}$ and $\tilde{Z}''$.

As shown in Figure 2, the computational mesh utilized here is more refined in the region of the domain close to the fuel jet, as this represents the region of interest. To reduce computational costs, both in the co-flow region and in the farthest zone from the flame, the mesh elements are bigger. The computational mesh employed features 33,200 elements. The smallest mesh element has a characteristic size of 0.125 mm and is located at the inner wall of the fuel discharge tube (exit region of the fuel jet), ensuring a $y+$ value smaller than 4 on the wall. Conversely, the largest mesh element, featuring a characteristic size of 4 mm, is located at the farthest radial position from the flame.
<figure><img src="image_3.png" /><figcaption>Figure 2 Computational mesh generated using the blockMesh utility within OpenFOAM, front and top view.</figcaption></figure>

To verify the obtained results, a thorough mesh independence study was performed, employing successively refined grids featuring 12,200, 22,300, and 33,200 elements. This study scrutinized temperature axial profiles and mass fractions of primary soot precursor species for each mesh configuration. Detailed findings from the referred mesh independence study are described in [52]. Overall, the associated analysis performed revealed that employing 22,300 mesh elements or more yielded negligible alterations in the results. Thus, to ensure consistency and reliability of the results, all simulations carried out in this work utilized a computational mesh comprising 33,200 elements.

**6. Results and Discussions**

The main results obtained in this work are outlined and examined in this section. More specifically, the flamelet libraries employed here are firstly described. This is followed by discussions of mean velocity and velocity fluctuations profiles. Next mixture fraction
and flame temperature results are presented and discussed as well. In the final part of
this section, PAH profiles and soot precursors are analyzed, and a comparative
assessment of soot volume fraction profiles is conducted.

6.1. Flamelet libraries

Firstly, flamelet libraries were generated using the libOpenSMOKE tool [49]. The detailed ABF chemical kinetic mechanism [45] in CHEMKIN format was used for this purpose. This mechanism was chosen because it a is widely employed one in soot formation studies [47]. Notice however that, when predicting the growth of PAH molecules, the ABF mechanism presents some limitations as it does not consider some important chemical pathways for PAH growth [53]. The ABF chemical kinetic mechanism features 101 chemical species, 543 chemical reactions, and the largest PAH accounted for is Pyrene (A4). Each flamelet was solved in this work using 82 non-equidistant points of stoichiometric scalar dissipation rate, covering a range of values from 10<sup>-6</sup> s<sup>-1</sup> (near-equilibrium) to 340 s<sup>-1</sup> (extinction). For an accurate representation of each flamelet, a total of 100 points along the mixture fraction were considered.

Accordingly, Figure 3a illustrates the relationship between the stoichiometric adiabatic
flame temperature and the stoichiometric scalar dissipation rate. This figure highlights
in particular the stable burning branch and the completely extinguished state. Notice
here that the transition point between these two states occurs at a stoichiometric scalar
dissipation rate corresponding to the quenching limit, approximately equal to 340 s<sup>-1</sup>.

Figure 3b shows in turn adiabatic flamelet temperature profiles as a function of Z at different χ<sub>st</sub>. This figure highlights the fact that the flamelets extinguish upon reaching the quenching limit stretch rate, resulting in inert mixing of the reactants.
To illustrate the adiabatic flamelet libraries created, the acetylene mass fraction corresponding to a $\chi_{st}$ of 0.05 s⁻¹ is compared in Figure 3c against the one obtained in reference [24] for the same flame configuration and the same $\chi_{st}$. As observed from this last figure, the agreement between the acetylene mass fraction profiles is relatively good, even though different chemical kinetic mechanisms were employed to obtain these results. Indeed, the past result [24] used for comparison purposes here was obtained with the NBP chemical kinetic mechanism, whereas the one in this work came from the ABF chemical kinetic mechanism.

Notice as well that, to account for radiative effects, following reference [40], about 10 distinct enthalpy defects ranging from 0 to -500 kJ/kg were employed. These generated flamelets were integrated using presumed Probability Density Function (PDF) models for mixture fraction, enthalpy defect, and scalar dissipation rate, and subsequently utilized by the RANS flamelet solver employed here.

<figure><img src="image_4.png" /><figcaption>Figure 3 Adiabatic steady flamelet solutions. a) Temperature as a function of $\chi_{st}$ at stoichiometric mixture fraction. b) Temperature as a function of Z at different $\chi_{st}$. c) Acetylene mass fraction at $\chi_{st} = 0.05$ s⁻¹ compared with literature results [24].</figcaption></figure>

## 6.2. Mean velocity and velocity fluctuations near the fuel jet exit

The characteristics of the turbulent velocity fields at the vicinity of the fuel jet exit is first analyzed here. This analysis is carried out using the DSM model because, as it will
be shown in Section 6.6, this model leads to the best predictions in terms of soot volume
fraction profiles. Accordingly, Figure 4 shows radial profiles of mean streamwise
velocity and its root-mean-square (RMS) at a flame axial position of 2 mm above the
burner. For comparison purposes, this figure also includes the experimental data
characterizing the Adelaide ENH1 jet flame [23] and other past numerical results
obtained using an LES approach [24]. From Figure 4a it is particularly observed that,
regarding the peak value on the centerline, in comparison to the experimental data, the
streamwise velocity obtained in this work exhibits discrepancies of less than 5 m/s
(~5%). In addition, this mean velocity profile closely aligns with the outcomes from
the previous LES simulation [24] accounted for, as evidenced by the agreement in their
peak values. This concordance between the RANS results obtained here and the
previous LES findings is noteworthy. Notice however that the velocity profile predicted
in this work is somehow more elongated in the radial direction in comparison to both
the experimental data and the past LES results. This difference may be attributed to the
inherent limitations of the RANS approach employed here, as it seems that the shear
layer, between 0.4 and 0.8 r/D, is not properly described here. Indeed, in this work, a
velocity of about 0 m/s is reached at 0.8 r/D, whereas at the experimental and past LES
results this velocity value is reached at around 0.5 r/D.

Similarly, the root-mean-square (RMS) velocity profile obtained here (Figure 4b) agrees
with the experimental data relatively well. It is worth noticing however that there is a
displacement from the centerline by a factor of 0.2 r/D. This displacement agrees with
the findings of the past LES study [24] used as reference here, which suggest that the
turbulent flow within the Adelaide Flame 1 burner may not have fully developed during
the experiments. These discrepancies in velocity fields led the past LES simulations to
predict a lifted flame, in contrast to the attached one observed in this study. This last
aspect is illustrated in Figure 4c, which shows the radial temperature profile at 2 mm above the burner.

Furthermore, it is worth noticing as well that the peak velocity fluctuation value obtained in this work matched the experimental data within 5%. These results outperformed the past LES study [24] accounted for here, which overestimated this fluctuation by nearly two times the experimental value, possibly due to the lifted flame configuration obtained there. Despite the observed discrepancies in the numerical results obtained in this work when compared to the experimental data, it is important to emphasize that those near the burner nozzle are unlikely to exert a significant influence further downstream, where the presence of soot becomes more prominent.

<figure><img src="image_5.png" /><figcaption>Figure 4 Radial profiles of a) mean velocity, b) velocity RMS, and c) mean temperature at a flame axial position of 2 mm above the burner.</figcaption></figure>

### 6.3. Mean mixture fraction and flame temperature

In the left plot of Figure 5, the mean mixture fraction obtained in this work along the flame centerline is compared with previous numerical results [24]. Significantly, a marked agreement becomes evident in the downstream part of the flame. Nevertheless, some discrepancies are observed at the axial positions where $z/D < 50$, primarily due to the contrast between the results obtained in the LES simulations [24], characteristics
of a lifted flame, and the experimental flame configuration, which represents an
attached one. Consequently, the present mixture fraction results align more closely with
those of the experimental measurements, given that the actual flame exhibits an attached
flame configuration as the one predicted here (Figure 4c).

In the right plot of Figure 5 in turn, the centerline flame temperature computed here, with and without radiation effects, is compared with the corresponding experimental data. This particular result highlights the accuracy of the numerical model employed here to predict flame temperatures when radiation effects are included. As observed in this figure indeed, the model closely captures the temperature peak value of about 1800K. However, a noteworthy observation is that, in the numerical predictions carried out in this work, there is a slight downstream shift of the peak temperature, suggesting an elongated flame structure or radiation heat transfer to the flame surroundings that has not been accounted for.

Furthermore, Figure 5 (right plot) draws attention to two scenarios, one incorporating
radiation effects and the other one excluding them. This differentiation underscores the
importance of accounting for radiation in sooting flames. Specifically, the exclusion of
radiation effects leads to a significant increase in the peak temperature by nearly 300 K.
Furthermore, the entire flame temperature profile increases in size, bringing it into a
disagreement with the experimental data. As observed in this figure, the previous LES
study [24] predicts the flame temperature fairly accurately as well.
<figure><img src="image_6.png" /></figure>

Figure 5 Mean mixture fraction (left) and flame temperature (with and without radiation) (right) profiles along the flame centerline compared against experimental data and past LES results [24].

## 6.4. PAH precursor analysis

A critical aspect of the numerical simulations carried out in this work is the accurate prediction of the PAH, which serve as the primary source for soot. Indeed, within the near field of jet flames, inception and condensation play a dominant role in soot evolution. Therefore, it is imperative to properly reproduce the distribution of PAH in this region to carry out accurate soot predictions.

In this context, as an attempt to improve PAH mass fraction predictions, in many past works [6,9], additional transport equations have been introduced into flamelet models.

In this work however, when studying the detailed soot models accounted for here, the PAH fields are derived from the steady laminar flamelet model only. This is done even though this combustion model has been shown to encounters challenges. Indeed, as highlighted in the DNS work by Bisetti et al. [54], which shows that, due to the relatively slow PAH chemistry, the PAH mass fraction deviates significantly from those computed using a steady laminar flamelet model. Owing to this model shortcoming and the lack of soot precursors species experimental data for this flame configuration, for
the sake of simplicity, the smallest aromatic species, benzene has been selected here as
the soot precursor in all detailed soot models. Thus, aiming to compare the effect of this
choice with other potential soot precursor species, Figure 6 gives centerline profiles
obtained for different PAH, namely benzene ($C_6H_6$), naphthalene ($C_{10}H_8$), and pyrene
($C_{16}H_{10}$), and also the profile characterizing the lumped PAH considered as soot
precursor in a past LES work [24]. Accordingly, for $0 \le z/D \le 20$, pyrene exhibits mass
fraction peak values in the order of $10^{-3}$, whereas naphthalene presents peaks in the
order of $10^{-5}$. The benzene mass fraction values are four times greater than those
characterizing the precursors used in the previous LES work [24]. Furthermore, the
profile shape of benzene differs from the bell-shaped profile observed in the referred
LES work. Since benzene is the main driver for soot formation in this work, these
substantial differences in PAH could inevitably manifest themselves in the results of
soot volume fraction computed in this work. Notice here that the assessment of the
influence on the obtained results of different PAH used as soot precursors is beyond the
scope of this work.

<figure><img src="image_7.png" /><figcaption>Figure 6 PAH mass fractions along the flame centerline compared against PAH used as soot precursors in past LES simulations [24]. Left plot: Close view of the fuel jet discharge zone. Right plot: Full $z/D$ extent.</figcaption></figure>
## 6.5. Soot volume fraction source terms analysis

In this section, the soot volume fraction source terms $df_v/dt$ obtained using all detailed soot models studied here (MOMIC, HMOM, and DSM) are presented and discussed. More specifically, to assess the accuracy and reliability of these soot models, their corresponding results are compared against numerical results obtained in prior studies (for inception + condensation) [24]. The analysis will be carried out by following the soot formation chain, from inception, condensation, surface growth to oxidation. It is worth noticing that, since there are many studies in the literature reporting the use of different inception efficiencies to lower the inception rate [48], the soot volume fraction results obtained in this work with all detailed soot models required the fine-tuning of the sticking coefficient $\gamma$, which acts as an inception efficiency. Thus, as shown in Table 1, for the MOMIC, HMOM, and DSM models, $\gamma$ values of $25 \times 10^{-3}$ in Eq. 14, $160 \times 10^{-6}$ in Eq. 23, and $40 \times 10^{-6}$ in Eq. 31 were employed, respectively. This tunning was performed by adjusting the computed maximum soot volume fraction along the centerline to the measured value to within 10% of discrepancy.

### 6.5.1. Source terms due to inception and condensation

Figure 7 (left) presents the computed soot volume fraction source terms $df_v/dt$ due to inception along the flame centerline. In this figure it is seen that the $df_v/dt$ computed values due to inception obtained using the MOMIC method are ten times larger than those computed with the HMOM, and nearly twenty times larger than that characterizing the DSM. These $df_v/dt$ related discrepancies are primarily attributed to the significantly higher value of the sticking coefficient used in MOMIC ($\gamma = 25 \times 10^{-3}$), which is more than two orders of magnitude greater than those used in HMOM ($\gamma = 160 \times 10^{-6}$) and DSM ($\gamma = 40 \times 10^{-6}$) (Table 1). The observed discrepancies also arise from the distinct inception models employed by each soot prediction method. Indeed,
the MOMIC model uses Eqs. (14) and (15), a concept similar to that of the DSM model, Eq. (31), whereas the HMOM model (Eq. (23)) is based on a completely different concept involving a prior dimerization process (PAH collision), with the inception term derived from the collision of these dimers. Note, also, that a substantially higher sticking coefficient is required for the MOMIC model to reproduce the experiments, because this model formulation does not account for condensation, and thus this sticking coefficient value enables to achieve a $df_v/dt$ from PAH (inception + condensation) similar to the values characterizing the HMOM and DSM, as will be discussed in section 6.5.2.

The $df_v/dt$ due to condensation along the flame centerline is plotted in Figure 7 (right). In this case, to avoid overprediction of soot and to achieve values comparable to those of HMOM, the DSM condensation source term (Eq. (32)) was reduced by a factor of 0.01 (Table 1). The HMOM source term remained unaltered and adhered to the its original formulation [9]. This was done because its condensation source term (Eq. (24)) already includes a reduction, as it is formulated based on the dimer concentration. Therefore, the DSM includes a reduction coefficient similar to the employed in the inception for HMOM. In Figure 7 (right), it can be observed that, in the case of HMOM, the source term due to condensation exhibits a peak value position similar to that of inception (Figure 7, left), around z/D=65. In contrast, in the DSM, condensation peaks at z/D=80, whereas inception does at z/D=62. This DSM peak axial position discrepancy results in a shift in the $df_v/dt$ due to PAH (condensation + inception), as shown in Figure 8, leading to a downstream-shifted soot volume fraction profile, which is analyzed in Section 6.6.
<figure><img src="image_8.png" /><figcaption>Figure 7 Soot volume fraction source terms, $df_V/dt$, along the flame centerline due to inception (left) and condensation (right).</figcaption></figure>

### 6.5.2. Source terms due to inception and condensation analysis

In Figure 8, $df_V/dt$ evolution along the flame centerline due to inception + condensation are shown, and compared with previous numerical results [24]. Notice that, unlike the present work that uses benzene obtained from the flamelet model as soot precursor, the previous numerical study [24] accounts for, as soot inception precursor, a group of PAH obtained by solving an additional transport equation. In addition, for the sake of comparison, the past numerical data [24] shown in Figure 8 is multiplied by 10. Figure 8 shows, in particular, that although the MOMIC model does not include condensation in its formulation, it achieves a peak value similar to that of the DSM, whose peak value is around 6.8 ppm/s. This similarity is primarily because the DSM features a very low value of $df_V/dt$ due to inception (peak of 0.2 ppm/s), and a relatively high value of $df_V/dt$ due to condensation (peak of 6.5 ppm/s), similar to the $df_V/dt$ due to inception in MOMIC (peak of 6.2 ppm/s). In the case of HMOM, the inception + condensation source terms feature a profile including two peak values, one at around $z/D=50$ and the other at around $z/D=65$ (Figure 8). This source term profile is similar to the benzene mass fraction one (Figure 6), underscoring that the inception +
condensation source term in this model follows the shape of the employed precursor profile. The inception source terms in MOMIC and DSM also show a profile featuring two peaks (Figure 7), although with different values which attributed to the sticking coefficient values employed.

In the case of DSM, the combined $d_f_v/dt$ due to inception + condensation does not exhibit two peak values. Instead, it shows a profile identical to that of the $d_f_v/dt$ due to condensation (Figure 7 right), this is mainly due to the relatively low values of $d_f_v/dt$ due to inception when compared to condensation. Some of the discrepancies between the $d_f_v/dt$ values computed here and those from the referred previous study arise from the fact that the peak values of the detailed soot models studied here have been adjusted using the sticking coefficient to align with soot volume fraction experimental values [23]. The differences between the HMOM and DSM profiles can be attributed to the different formulations employed in each model. For instance, in DSM, there is an inter-sectional soot dynamics allowing soot to pass from one section to another due to condensation. While DSM and MOMIC display similar peak values, HMOM shows a peak value of 3 ppm/s, which is half that characterizing the other two models. This is attributed to the fact that HMOM presents a significantly higher value of $d_f_v/dt$ due to surface growth compared to the other cases, which is given in Figure 9 (left) and that will be further discussed in section 6.5.3.

Regarding the $d_f_v/dt$ due to inception + condensation obtained previously in reference [24], although its peak value is about ten times smaller to those predicted in this work, it exhibits a profile similar to that of DSM, possibly because its lumped PAH profile does not show the two peaks observed in this study (Figure 6). Compared to HMOM and MOMIC ($z/D=65$), its peak value is however slightly downstream, at $z/D=70$, whereas the DSM peak is further downstream (at $z/D=75$), due to its predominance of the
condensation source term compared to inception. It should be noted that the numerical
work [24] used for comparison purposes employed the HMOM for soot modelling also,
which may explain why the peak value of the HMOM source term is located at the
lumped PAH peak around z/D=70 (Figure 8). In contrast, the DSM model predicts a
profile shifted downstream from the peak value of the lumped PAHs utilized in [24].

<figure><img src="image_9.png" /><figcaption>Figure 8 Soot volume fraction source terms $df_v/dt$ along the flame centerline due inception + condensation.</figcaption></figure>

6.5.3. Source terms due to surface growth and oxidation analysis

Profiles of $df_v/dt$ due to surface growth and oxidation along the flame centerline are shown in turn in Figure 9. Regarding surface growth first, although they present different peak values, HMOM 3.9 ppm/s and MOMIC 1.3 ppm/s, the two moment-based models employed here feature a similarly positioned bell-shaped profile. This can be attributed to the fact that these terms do not depend on the precursor PAH but on acetylene. Additionally, the surface growth term in MOMIC has not been adjusted and it retains its original formulation, whereas in the case of HMOM, the rate constant $k_{sg}$ in Eq. (25) has been adjusted by a factor of 0.002 (Table 1). In addition, in the DSM, to reduce the excessive soot increase initially caused by surface growth, a factor of 0.014
was included in $k_4$, Eq. (33) (Table 1). In other works, the DSM surface growth has also been modified to match experimental data [45,47]. These adjustments were made to compare these source term values with those characterizing the MOMIC model.

For DSM, the surface growth related profile (Figure 9, left) is shifted downstream compared to the other two models. This trend was also observed in the $df_v/dt$ profile due to condensation (Figure 7, right), which can be attributed to the way the source terms are modeled in DSM. Indeed, as highlighted in Section 6.5.2, the DSM method includes an inter-sectional soot dynamics approach that the moment-based methods do not. This could be one of the reasons why both the condensation and surface growth terms are shifted downstream. Regarding the numerical values obtained using each detailed soot model, HMOM presents the largest values of $df_v/dt$ due to surface growth (peak of 3.8 ppm/s). Compared to the 3 ppm/s peak observed in the $df_v/dt$ due to condensation + inception, it can be argued that the main source of soot in this model is due to surface growth (accounting for the adjustments made to the models). The differences in the peaks values of $df_v/dt$ due to inception + condensation observed in Figure 8 between the DSM and MOMIC models and the HMOM one, nearly doubling the first models to the latter one, is somehow compensated by the relatively high values of surface growth associated with the HMOM (Figure 9, left). Indeed, in this figure the MOMIC model yields a peak value of 1.3 ppm/s, and the DSM, a maximum value of 0.9 ppm/s, which is small compared to the $df_v/dt$ due to condensation + inception (peak of 6.8 ppm/s).
<figure><img src="image_10.png" /></figure>

Figure 9 Soot volume fraction source terms $df_v/dt$ along the flame centerline due to surface growth (left) and oxidation (right).

Concerning the $df_v/dt$ due to oxidation along the flame centerline (Figure 9, right), the three soot models studied here show quite different profiles. In the case of DSM, oxidation begins to increase after $z/D=100$, whereas in the other two models, oxidation starts at $z/D=70$. Consequently, the DSM model is expected to lead a more elongated soot profile, Section 6.6. As highlighted in Table 1, in the HMOM, the oxidation term $k_{ox}$ present in Eq. 26 was adjusted by a factor of 0.1, while in the DSM $k_5$ appearing in Eq. 34 was reduced by a factor of 0.02, and no adjustments were made in the MOMIC model, so the associated oxidation term remained the same as in its original formulation. These adjustments were made because the original excessive oxidation values reduce the extent of the final soot profile. Additionally, these adjustments were necessary to obtain oxidation rate values comparable with those obtained with the MOMIC model. Accordingly, the MOMIC model presents the smallest oxidation source term absolute value, with an absolute value of 0.05 ppm/s. The HMOM model features in turn the highest, with a peak value of 0.7 ppm/s. The influence of these oxidation terms will now be seen in the radial profiles of the $df_v/dt$ global source terms.
6.5.4. Global (inception, condensation, surface growth and oxidation) source terms analysis

In Figure 10 the global source terms, df<sub>v</sub>/dt, including inception, condensation, surface growth, and oxidation, are shown along the flame centerline and at three different axial positions, z/D=98, z/D=103, and z/D=117. It is worth emphasizing that the results shown in this figure do not include soot growth due to coagulation, as coagulation does not change the soot volume fraction. From this figure, it is observed first that, along the centerline, all three studied models exhibit a similar bell shape and a similar peak value of about 7 ppm/s. Additionally, in the MOMIC and HMOM models, the peak value occurs between 60≤z/D≤70. In contrast, due to the downstream shift in condensation and surface growth in this model (which are the primary contributors to the global source term), the DSM peaks at z/D=80. The HMOM profile occurs slightly downstream of the MOMIC one because surface growth has a greater influence than inception. In MOMIC, the primary contributor to soot growth is inception, which profile is shape similar to that of the benzene precursor.

In the MOMIC model, oxidation has a minimal effect on the global source term along
the axial direction, whereas in HMOM, negative values of the global source term appear
between 110≤z/D≤140. In DSM, negative values appear downstream z/D=140 only.
This last outcome will be shown to significantly affect the computed soot volume
fraction distribution, making it more elongated in the DSM model. As seen in Figure 10,
the influence of oxidation on the obtained results is more evident along the flame radial
direction. For instance, at z/D=76, near the peak soot concentration for all three models,
the MOMIC model shows significant oxidation rates with values close to -30 ppm/s. In
contrast, the other two soot models show minimal oxidation rates, with HMOM
showing values around -1 ppm/s and DSM featuring no oxidation at all. These
discrepancies in the obtained results are due to the adjustment factors used in the
oxidation terms for DSM and HMOM, the inherent formulation of the soot models, and
the fact that oxidation begins further downstream in DSM (z/D=140). At the other two
axial positions analyzed here, z/D=98 and z/D=103, the trends are similar. Indeed, the
MOMIC model shows values greater than -5 ppm/s, whereas HMOM features peak
values around -1 ppm/s. The DSM model does not show oxidation rates because, as
seen in the top left plot of Figure 10, oxidation in this case begins at z/D=140.
Accordingly, this model does not exhibit significant oxidation rates at the axial
positions accounted for. The oxidation effect along the flame radial direction narrows
the soot profiles. Consequently, as will be shown now, the narrowest soot volume
fraction profile is observed for MOMIC, followed by HMOM, with the widest soot
profile being the DSM one.
<figure><img src="image_11.png" /><figcaption>Figure 10 Soot volume fraction global source terms, $df_v/dt$, along the flame centerline and along the radial direction at different axial positions, $z/D = 76, 98, \text{and } 103$.</figcaption></figure>

## 6.6. Soot volume fraction profiles analysis

Profiles of soot volume fraction (SVF) along the flame centerline are shown in Figure 11. For comparison purposes, this figure also includes the corresponding experimental data [23] and past LES results [24]. Due to the tuning of the sticking factors analyzed in Section 6.5, the SVF predicted peak values agree well with the experimental ones. However, there are some discrepancies in the predicted axial locations of these peak values for the 2EQ, MOMIC, and HMOM models. The locations predicted by these soot models lie significantly upstream, at around $z/D = 85$, whereas the corresponding experimental data indicate $z/D = 105$. The model that best predicts the position of the soot peak value is the DSM, at $z/D = 100$.
Concerning the detailed soot models, Figure 11 indicates that the initial increase in soot volume fraction occurs first for the MOMIC, then for the HMOM, and finally for the DSM. Since the primary source of soot growth in MOMIC is inception, the source terms analyses carried out in previous sections explains this result obtained for MOMIC, which presents between $15 \le z/D \le 35$ a global source term larger than the other two detailed soot models (Figure 10). These values could be of course adjusted by reducing the influence of inception and increasing growth due to surface growth, which would require another parameter adjustment. The next model to show an increase in soot volume fraction is the HMOM, and the main source of growth in this model is condensation, and then surface growth. Finally, the DSM is the model that best captures the initial phase of soot formation, since both condensation (Figure 7) and surface growth (Figure 9) source terms profiles are shifted downstream compared to the MOMIC and HMOM ones. Downstream the peak values, the soot volume fraction decreases along the centerline at a similar rate for HMOM and MOMIC and at a slower rate for the DSM. This outcome is primarily attributed to the fact that the oxidation rate in DSM occurs much later compared to the other models, i.e., downstream $z/D \approx 140$. Consequently, an adjustment or change in the oxidation rate would be necessary to better predict the final soot profile with the DSM model, one that has not been attempted here.

It is also interesting to notice the results obtained with the 2EQ model, which uses acetylene as its soot precursor. The 2EQ model utilized in this work remains consistent with its original formulation, with all constants unaltered from the specifications of the model. Remarkably, even without any empirical constant adjustment, this 2EQ model offers relatively accurate predictions of the peak soot value, with a discrepancy of only 75 ppb. However, it is essential to emphasize that, qualitatively, this model yields the
smallest soot profile axial extent, i.e., 40 z/D, which differs from the experimental data, where the soot profile has an axial extent of about 100 z/D. Thus, this 2EQ model deviates the most from the experimental measurements.

Conversely, for the soot detailed models that use benzene as the soot precursor species, the centerline soot profiles exhibit a remarkable similarity for all models, as also observed in the global source term profiles obtained (Figure 10, top left). The soot regions have a similar axial extent except for the DSM, which has lower oxidation rate in the final part of the flame, resulting in a more elongated soot profile. Among all the detailed soot models studied here, in terms of soot volume fraction, the DSM shows the closest agreement with the experimental data. In addition, the DSM is the only model that captures the position of the experimental soot peak value, which can be attributed to the shifted condensation (Figure 7) and surface growth (Figure 9) source term profiles compared to the moment-based models studied here.

It is also worth noticing that several discrepancies in the soot profiles obtained with the detailed soot models accounted for are evident in the final part of the domain, where the oxidation process model plays a pivotal role in either elongating or reducing the final soot volume fraction profile. The DSM model, characterized by its conspicuously elongated soot profile, matches the overall soot volume fraction profile obtained in the experiments [23], albeit with an elongated tail. It is anticipated that by fine-tuning and increasing the oxidation rate, the final soot volume fraction profile could become less elongated, thus approaching experimental values more closely. Therefore, further adjusting the empirical constants within these soot models has the potential to yield results that better match the experimental observations.
<figure><img src="image_12.png" /></figure>

Figure 11 Soot volume fraction profiles along the flame centerline compared with experimental data [23] and past LES results [24].

Finally, radial profiles of soot volume fraction at different flame axial positions are shown in Figure 12. From this figure, it is observed that the DSM model provides the best description of the radially distributed soot. This is primarily because the DSM model features an elongated soot volume fraction axial profile (Figure 11), and the radial positions assessed in Figure 12 extend from the peak of soot volume fraction to the final part of the flame ($z/D > 89$). It should be noticed that the radial amplitude of the soot profiles is correctly predicted by the DSM model in the first three axial positions compared here. However, it deviates significantly, nearly doubling the experimental data in later positions. This deviation is mainly due to the smaller oxidation source term absolute value in this model. The MOMIC model deviates the most from the experimental data, possibly due to a higher radial oxidation rate (Figure 10), and to the fact that the peak soot value in this method is located at around $z/D \approx 80$, whereas the plots in Figure 12 are for $z/D > 89$. The HMOM model is in turn the second-best model in capturing the radial profiles. This is because the soot peak value
of the HMOM is located slightly downstream of the MOMIC one, originated from the fact that, unlike MOMIC that primarily grows soot through inception, HMOM presents increased surface growth rates. Also, HMOM exhibits axial oxidation rates smaller than MOMIC.

<figure><img src="image_13.png" /><figcaption>Figure 12 Radial soot volume fraction profiles at different flame axial positions, compared with experimental data [23].</figcaption></figure>

## 7. Conclusions

In this work, an assessment of four distinct soot models (namely, 2EQ, MOMIC, HMOM, and DSM) was conducted. The referred assessment was carried out within the context of RANS simulations of turbulent non premixed flames. The outcomes from this evaluation were systematically compared against both experimental data and past LES results available in the literature. As a part of the assessments undertaken here, an
initial comparison of mean velocity and velocity fluctuations close to the burner feed tube was performed. The associated findings indicate that the RANS model, augmented with the modifications applied to the standard k-ε model, successfully approximated the velocity profile in close agreement with the LES results. However, the velocity peak value fell below the experimental measurements by about 5 m/s (~5%). Regarding the velocity fluctuations, the model succeeded in capturing the peak values observed in the experiments, although these peak values were located at a radial position r/D = 0.1 away from that characterizing the experimental data. The velocity fluctuations related results surpassed the performance of the LES work used as reference here, mainly because the numerical simulations performed here lead to an attached flame.

In terms of flame temperature, the numerical model employed here reproduced the experimental temperature profiles with a discrepancy of 10% only in the flame downstream region. In particular, the incorporation of radiation effects into the numerical model enabled to improve the representation of the temperature field. Regarding the selection of the PAH to be used as soot precursor species in the detailed soot models employed here, the decision was made to use benzene as soot precursor due to inherent limitations of the flamelet model in generating high-quality profiles for larger PAH species. This contrasts with the LES work used as reference in this work, which considered multiple PAH. It is important to notice that the benzene mass fraction profile used here proved to be suitable for accurately capturing the location of the peak PAH concentration, resembling in this way the results obtained in the LES study. Nonetheless, the benzene concentration values were approximately four times higher than the PAH ones used in the LES work.

The analysis of soot volume fraction source terms reveals discrepancies between the models. The MOMIC model showed inception source terms significantly higher than
the HMOM and DSM ones, attributed to its higher sticking coefficient and different inception models. The condensation source terms were discrepant, with DSM showing a further peak, when compared to HMOM. These differences led to variations in the combined inception and condensation source term profiles, impacting the soot volume fraction results. Adjustments in the sticking coefficients were used to align the numerical results obtained using these models with experimental data, highlighting the complexity of accurately modeling soot formation in turbulent flames.

When comparing the soot volume fraction profiles, it was found that the initial increase in soot volume fraction occurs first for MOMIC, then for HMOM, and finally for DSM. The dependence of MOMIC on inception as the primary source of soot growth led to an earlier increase, with HMOM and DSM showing an increased influence from condensation and surface growth. Downstream, the soot volume fraction decreased at different rates, with DSM showing a slower decrease due to its later oxidation onset. The 2EQ model, using acetylene as soot precursor, yielded the smallest soot profile and deviated the most from experimental data.

The soot volume fraction profiles along the flame centerline highlight that, despite agreeing on peak values, there are significant discrepancies in the predicted axial locations of these peaks. The DSM model best captured the axial location of experimental soot peak, whereas the other models predicted it upstream. Adjustments in the oxidation rate are suggested to improve the final soot profile predictions. Radial profiles of soot volume fraction further showed that the DSM model provided the best description of radially distributed soot, closely matching experimental data in the first three axial positions accounted for here, but deviating at later positions due to lower oxidation rates.
Overall, this work underscores the importance of selecting appropriate PAH precursors and fine-tuning model parameters to improve the accuracy of soot predictions. The DSM model emerged as the most promising one for capturing both axial and radial soot volume fraction profiles, though adjustments in oxidation rates seem to be necessary to better match experimental measurements.

## Acknowledgements

This work has been supported by CONCYTEC-FONDECYT (PROCIENCIA) (Peru), Contract No. 415-2019-2019-FONDECYT. Support was also received by the Brazilian Conselho Nacional de Desenvolvimento Cientifico e Tecnologico, CNPq, under Research Grants No. 403904/2016-1 and 304444/2018-9. For the purpose of Open Access, a CC-BY public copyright license has been applied by the authors to the present document and will be applied to all subsequent versions up to the Author Accepted Manuscript arising from this submission.

## Disclosure statement

No potential conflict of interest was reported by the authors.

## References

[1] J.W. Martin, M. Salamanca and M. Kraft, *Soot inception: Carbonaceous nanoparticle formation in flames*, Prog Energy Combust Sci 88 (2022), pp. 100956.

[2] D.D. Li, C. Wang, Q.N. Chan and G.H. Yeoh, *Soot: A review of computational models at different length scales*, Exp Comput Multiph Flow 5 (2023), pp. 1-14.
[3] S. Valencia, S. Ruiz, J. Manrique, C. Celis and L.F. Figueira da Silva, *Soot modeling in turbulent diffusion flames: review and prospects*, Journal of the Brazilian Society of Mechanical Sciences and Engineering 43 (2021), pp. 219.

[4] N. Peters, *Turbulent Combustion*, Cambridge University Press, 2000.

[5] M.E. Mueller and H. Pitsch, *LES model for sooting turbulent nonpremixed flames*, Combust Flame 159 (2012), pp. 2166-2180.

[6] S. Yang, J.K. Lew and M.E. Mueller, *Large Eddy Simulation of soot evolution in turbulent reacting flows: Presumed subfilter PDF model for soot-turbulence-chemistry interactions*, Combust Flame 209 (2019), pp. 200-213.

[7] H. Maldonado Colmán, A. Cuoci, N. Darabiha and B. Fiorina, *A virtual chemistry model for soot prediction in flames including radiative heat transfer*, Combust Flame 238 (2022), pp. 111879.

[8] F. Nmira, A. Bouffard and J.-L. Consalvi, *Flamelet/transported PDF simulations of ethylene/air jet turbulent non-premixed flame using a three-equation PAH-based soot production model*, Combustion Theory and Modelling 27 (2023), pp. 820-851.

[9] M. Mueller, G. Blanquart and H. Pitsch, *Hybrid Method of Moments for modeling soot formation and growth*, Combust Flame 156 (2009), pp. 1143-1155.

[10] S. Yang and M.E. Mueller, *A Multi-Moment Sectional Method (MMSM) for tracking the soot Number Density Function*, Proceedings of the Combustion Institute 37 (2019), pp. 1041-1048.

[11] I.M. Kennedy, *Models of soot formation and oxidation*, Prog Energy Combust Sci 23 (1997), pp. 95-132.
[12] M. Frenklach, *Method of moments with interpolative closure*, Chem Eng Sci 57 (2002), pp. 2229-2239.

[13] S. Salenbauch, A. Cuoci, A. Frassoldati, C. Saggese, T. Faravelli and C. Hasse, *Modeling soot formation in premixed flames using an Extended Conditional Quadrature Method of Moments*, Combust Flame 162 (2015), pp. 2529-2543.

[14] C.A. Hoerlle and F.M. Pereira, *Effects of CO2 addition on soot formation of ethylene non-premixed flames under oxygen enriched atmospheres*, Combust Flame 203 (2019), pp. 407-423.

[15] F. Gelbard, Y. Tambour and J.H. Seinfeld, *Sectional representations for simulating aerosol dynamics*, J Colloid Interface Sci 76 (1980), pp. 541-556.

[16] K. Netzell, H. Lehtiniemi and F. Mauss, *Calculating the soot particle size distribution function in turbulent diffusion flames using a sectional method*, Proceedings of the Combustion Institute 31 (2007), pp. 667-674.

[17] P. Rodrigues, B. Franzelli, R. Vicquelin, O. Gicquel and N. Darabiha, *Coupling an LES approach and a soot sectional model for the study of sooting turbulent non-premixed flames*, Combust Flame 190 (2018), pp. 477-499.

[18] M. Balthasar and M. Kraft, *A stochastic approach to calculate the particle size distribution function of soot particles in laminar premixed flames*, Combust Flame 133 (2003), pp. 289-298.

[19] H. Zhao, C. Zheng and M. Xu, *Multi-Monte Carlo method for coagulation and condensation/evaporation in dispersed systems*, J Colloid Interface Sci 286 (2005), pp. 195-208.

[20] A. Kalbhor, D. Mira and J. van Oijen, *A computationally efficient approach for soot modeling with discrete sectional method and FGM chemistry*, Combust Flame 255 (2023), pp. 112868.
[21] S.J. Brookes and J.B. Moss, *Measurements of soot production and thermal radiation from confined turbulent jet diffusion flames of methane*, Combustion Flame 116 (1999), pp. 49–61.

[22] D.L. Marchisio and R.O. Fox, *Solution of population balance equations using the direct quadrature method of moments*, J Aerosol Sci 36 (2005), pp. 43–73.

[23] S. Mahmoud, G.J. Nathan, Z.T. Alwahabi, Z.W. Sun, P.R. Medwell and B.B. Dally, *The effect of exit Reynolds number on soot volume fraction in turbulent non-premixed jet flames*, Combust Flame 187 (2018), pp. 42–51.

[24] S. Yang, J.K. Lew and M.E. Mueller, *Large Eddy Simulation of soot evolution in turbulent reacting flows: Strain-Sensitive Transport Approach for Polycyclic Aromatic Hydrocarbons*, Combust Flame 220 (2020), pp. 219–234.

[25] T. Zhang, L. Zhao, M.R. Kholghy, S. Thion and M.J. Thomson, *Detailed investigation of soot formation from jet fuel in a diffusion flame with comprehensive and hybrid chemical mechanisms*, Proceedings of the Combustion Institute 37 (2019), pp. 2037–2045.

[26] M. Frenklach and H. Wang, *Detailed modeling of soot particle nucleation and growth*, Symposium (International) on Combustion 23 (1991), pp. 1559–1566.

[27] X. Mercier, O. Carrivain, C. Irimiea, A. Faccinetto and E. Therssen, *Dimers of polycyclic aromatic hydrocarbons: the missing pieces in the soot formation process*, Physical Chemistry Chemical Physics 21 (2019), pp. 8282–8294.

[28] M. Frenklach and A.M. Mebel, *On the mechanism of soot nucleation*, Physical Chemistry Chemical Physics 22 (2020), pp. 5314–5331.

[29] G. Blanquart and H. Pitsch, *Analyzing the effects of temperature on soot formation with a joint volume-surface-hydrogen model*, Combust Flame 156 (2009), pp. 1614–1626.
[30] M. Mueller, G. Blanquart and H. Pitsch, *A joint volume-surface model of soot aggregation with the method of moments*, Proceedings of the Combustion Institute 32 I (2009), pp. 785–792.

[31] K.M. Pang, N. Karvounis, J.H. Walther, J. Schramm, P. Glarborg and S. Mayer, *Modelling of temporal and spatial evolution of sulphur oxides and sulphuric acid under large, two-stroke marine engine-like conditions using integrated CFD-chemical kinetics*, Appl Energy 193 (2017), pp. 60–73.

[32] S.H. Park and S.N. Rogak, *A novel fixed-sectional model for the formation and growth of aerosol agglomerates*, J Aerosol Sci 35 (2004), pp. 1385–1404.

[33] A. Bouaniche, L. Vervisch and P. Domingo, *A hybrid stochastic/fixed-sectional method for solving the population balance equation*, Chem Eng Sci 209 (2019), pp. 115198.

[34] B.E. Launder and B.I. Sharma, *Application of the energy-dissipation model of turbulence to the calculation of flow near a spinning disc*, Letters in Heat and Mass Transfer 1 (1974), pp. 131–137.

[35] B.B. Dally, A.R. Masri, R.S. Barlow and G.J. Fiechtner, *Instantaneous and Mean Compositional Structure of Bluff-Body Stabilized Nonpremixed Flames*, Combust Flame 114 (1998), pp. 119–148.

[36] A. Rowhani, Z.W. Sun, P.R. Medwell, Z.T. Alwahabi, G.J. Nathan and B.B. Dally, *Effects of the Bluff-Body Diameter on the Flow-Field Characteristics of Non-Premixed Turbulent Highly-Sooting Flames*, Combustion Science and Technology 194 (2022), pp. 378–396.

[37] N. Peters, *Laminar diffusion flamelet models in non-premixed turbulent combustion*, Prog Energy Combust Sci 10 (1984), pp. 319–339.

[38] K. Bray and N. Peters, *in*, *Turbulent Reacting Flows (1994),*.
[39] F. Liu, K.A. Thomson, H. Guo and G.J. Smallwood, *Numerical and experimental study of an axisymmetric coflow laminar methane-air diffusion flame at pressures between 5 and 40 atmospheres*, Combust Flame 146 (2006), pp. 456-471.

[40] A. Cuoci, A. Frassoldati, T. Faravelli and E. Ranzi, *Kinetic modeling of soot formation in turbulent nonpremixed flames*, Environ Eng Sci 25 (2008), pp. 1407-1422.

[41] R.S. Barlow, A.N. Karpetis, J.H. Frank and J.Y. Chen, *Scalar profiles and NO formation in laminar opposed-flow partially premixed methane/air flames*, Combust Flame 127 (2001), pp. 2102-2118.

[42] M.V. Smoluchowski, *Mathematical theory of the kinetics of coagulation of colloidal systems*, Z. Phys. Chem. 92 (1917), pp. 129-168.

[43] M. Balthasar and M. Frenklach, *Detailed kinetic modeling of soot aggregate formation in laminar premixed flames*, Combust Flame 140 (2005), pp. 130-145.

[44] A. Wick, M. Frenklach and H. Pitsch, *Systematic assessment of the Method of Moments with Interpolative Closure and guidelines for its application to soot particle dynamics in laminar and turbulent flames*, Combust Flame (2020), .

[45] J. Appel, H. Bockhorn and M. Frenklach, *Kinetic modeling of soot formation with detailed chemistry and physics: laminar premixed flames of C2 hydrocarbons*, Combust Flame 121 (2000), pp. 122-136.

[46] F. Gelbard and J.H. Seinfeld, *Numerical solution of the dynamic equation for particulate systems*, J Comput Phys 28 (1978), pp. 357-375.

[47] M.J. Thomson, *Modeling soot formation in flames and reactors: Recent progress and current challenges*, Proceedings of the Combustion Institute 39 (2023), pp. 805-823.
[48] M. Saffaripour, A. Veshkini, M. Kholghy and M.J. Thomson, *Experimental investigation and detailed modeling of soot aggregate formation and size distribution in laminar coflow diffusion flames of Jet A-1, a synthetic kerosene, and n-decane*, Combust Flame 161 (2014), pp. 848-863.

[49] A. Cuoci, A. Frassoldati, T. Faravelli and E. Ranzi, *OpenSMOKE++: An object-oriented framework for the numerical modeling of reactive systems with detailed kinetic mechanisms*, Comput Phys Commun 192 (2015), pp. 237-264.

[50] S.V. Patankar and D.B. Spalding, *A calculation procedure for heat, mass and momentum transfer in three-dimensional parabolic flows*, Int J Heat Mass Transf 15 (1972), pp. 1787-1806.

[51] N.H. Mohamed Ibrahim and M. Udayakumar, *Coupling RSM with soot model for the study of soot formation in a momentum-dominated strained jet flames*, J Therm Anal Calorim 141 (2020), pp. 2369-2389.

[52] S. Valencia, C. Celis and L.F. Figueira da Silva, *Application of a soot formation model based on an interpolative closure Method of Moments to a turbulent non-premixed flame*, in *Proceedings of the 26th International Congress of Mechanical Engineering*, 2021.

[53] A. Khrabry, I.D. Kaganovich, Y. Barsukov, S. Raman, E. Turkoz and D. Graves, *Compact and accurate chemical mechanism for methane pyrolysis with PAH growth*, Int J Hydrogen Energy 56 (2024), pp. 1340-1360.

[54] F. Bisetti, G. Blanquart, M.E. Mueller and H. Pitsch, *On the formation and early evolution of soot in turbulent nonpremixed flames*, Combust Flame 159 (2012), pp. 317-335.