A Detailed PAH and Soot Model for Complex Fuels in CFD Applications

Abstract

A model to predict soot evolution during the combustion of complex fuels is presented. On one hand, gas phase, $\hbox {polycyclic aromatic hydrocarbon (PAH)}$ and soot chemistry are kept large enough to cover all relevant processes in aero engines. On the other hand, the mechanisms are reduced as far as possible, to enable complex computational fluid dynamics (CFD) combustion simulations. This is important because all species transport equations are solved directly in the $\hbox {CFD}$. Moreover, emphasis is placed on the applicability of the model for a variety of fuels and operating conditions without adjusting it. A kinetic scheme is derived to describe the chemical breakdown of short- and long-chain hydrocarbon fuels and even blends of them. $\hbox {PAHs}$ are the primary soot precursors which are modeled by a sectional approach. The reversibility of the interaction between different $\hbox {PAH}$ classes is achieved by the introduction of $\hbox {PAH}$ radicals. Soot particles are captured by a detailed sectional approach too, which takes a non-spherical growth of particles into account. In this way the modeling of surface processes is improved. The applicability and validity of the gas phase, $\hbox {PAH}$, and soot model is demonstrated by a large number of shock tube experiments, as well as in atmospheric laminar sooting flames. The presented model achieves excellent results for a wide range of operating conditions and fuels. One set of model constants is used for all simulations and no case-dependent optimization is required.

1 Introduction

Low emissions during combustion of hydrocarbon-based fuels are a key criterion for the approval of aero engines. Along with carbon and nitrogen oxides, soot pollution is subject to research in the field of environment, climate and health (Ammann et al. 1998; Prüss-Ustün et al. 2016). The occurrence of cirrus clouds due to condensation trails is significantly enhanced by soot emissions from aero engines. These are only some of many consequences arising from undesired combustion products (Burkhardt et al. 2018; Travis et al. 2002). From a technical perspective, soot radiation in aero engines may cause unwantedly high wall heat loads (Wulff and Hourmouziadis 1997). For these reasons and the upcoming more stringent legal regulations, combustion concepts with higher efficiencies and less pollutant emissions are essential. This requires reliable soot prediction models able to reproduce its complex dependencies on temperature, pressure, fuel, and level of premixing.

Early empirical soot models correlated the soot volume fraction with known state variables, for example the temperature or the amount of unburned fuel (Edelmann and Harsha 1978). Based on such models, semi-empirical variants have been derived which take more physical effects into account. As in any model there are a number of constants which are inherently required. A disadvantage of simple soot models is given by constants which are only valid for the application they have been tuned for. To avoid this constraint more physical and chemical processes of soot formation have to be included. $\hbox {PAHs}$ were identified to be the major soot precursors, and the $\hbox {hydrogen abstraction and carbon addition (HACA)}$ mechanism was developed (Frenklach et al. 1985). Today, this approach is the basis of most detailed soot models (Kazakov et al. 1995; Wang and Frenklach 1997; Bisetti et al. 2012; Dworkin et al. 2011). Three approaches for soot are found most frequently in literature: the method of moments, the sectional approach, and the two-equation model. The method of moments approximates the $\hbox {soot particle size distribution (PSD)}$ by a predefined number of moments (Frenklach and Harris 1987). This method is very efficient but introduces a closure problem in the moments’ source terms (Mueller et al. 2009a). By the definition of suitable moment closures, the aggregation of soot particles can also be taken into account (Mueller et al. 2009b; Salenbauch et al. 2019). Aggregation describes the formation of non-spherically shaped soot particles, whereby surface chemical effects are more accurately reproduced. The sectional approach discretizes the $\hbox {PSD}$ by a finite number of sections, which also allows aggregation to be considered (Pope and Howard 1997; Richter et al. 2005; Dworkin et al. 2011; Rodrigues et al. 2018; Tian et al. 2021). A sufficiently high number of soot sections is required to accurately cover the complete soot mass spectrum. This increases the computational effort. In contrast, two-equation models (Leung et al. 1991; Di Domenico et al. 2010; Franzelli et al. 2015) neglect a detailed representation of soot surface effects and are computationally less expensive. Here, soot is modeled as a set of monodisperse particles and no information on the $\hbox {PSD}$ is provided. However, by introduction of an additional transport equation (three-equation model) (Franzelli et al. 2019), information on the $\hbox {PSD}$ can be obtained.

For aero engine soot simulations usually combustion models are employed based on tabulated chemistry (flamelet, FGM) (Koo et al. 2019; Mueller and Pitsch 2012). This keeps the computational effort low and allows flexibility with respect to turbulence modeling and computational grid. Another approach that is followed in this paper, is the direct solution of $\hbox {finite-rate chemistry (FRC)}$ in the $\hbox {CFD}$ solver. Here, a more or less large number of species transport equations has to be solved (Smooke et al. 1999; Zamuner and Dupoirieux 2000; Gerlinger 2005; Blacha et al. 2011; Dworkin et al. 2011; Franzelli et al. 2016; Veshkini et al. 2016; Pejpichestakul et al. 2019). This is computationally more expensive but offers the advantage of an inherently consistent treatment of all chemical species ranging from gas phase over the $\hbox {PAHs}$ to soot. In this way element and mass conservation is satisfied over the complete model, and it may be used for combustors with premixed and non-premixed regions. Moreover, the mixing with burned recirculated gas (what is typical for aero engines) does not violate model assumptions, as for example in case of the flamelet approach. Cooled walls or radiation losses are easy to implement, too. On the other hand, a relatively large kinetic scheme is required for a kerosene surrogate. Therefore, the development of an optimized, reduced gas phase mechanism is of particular importance, what is a part of this paper. A corresponding 3D $\hbox {Reynolds-Averaged Navier-Stokes (RANS)}$ soot simulation of an aero engine is demonstrated in Steinbach et al. (2018) and $\hbox {Large-eddy simulation (LES)}$ probably may be performed, too, if the mechanism provided in this paper is further reduced. This can be achieved by a test case depending reduction, e.g., for a constant combustor pressure.

The transport equations for soot can be solved by an Eulerian approach in analogy to gas phase species and $\hbox {PAHs}$, or using a Lagrangian tracking. The second method is often used in combination with flamelet models (Pitsch et al. 2000; Dellinger et al. 2020; Gallen et al. 2019).

To reduce the computational effort associated with the $\hbox {PAHs}$, a lumping of $\hbox {PAH}$ molecules into virtual sections is used (Frenklach and Wang 1991; Richter et al. 2005). Contrary to this, more detailed approaches consider each individual $\hbox {PAH}$ as an own species (Blanquart et al. 2009; Djokic et al. 2014). This causes a significantly higher computational effort. For example, in comparison with our model the mechanism of Blanquart et al. (2009) with 148 species requires $60\%$ more computational time in case of the laminar kerosene flame in Sect. 3.5 without soot model. To avoid such an increase and keep the computational effort low, lumping is reasonable for the $\hbox {PAHs}$ with the associated modeling uncertainties.

The majority of the current soot models is designed for simple hydrocarbons such as methane or ethylene (Carbonell et al. 2009; Eberle et al. 2015; Saggese et al. 2015). Complex fuels which are used for example in aero engines, cause an additional complexity in the prediction of soot. Highly detailed kinetic schemes have been developed for large hydrocarbon fuels (Park and Rogak 2004; Richter et al. 2005; Ciajolo et al. 2009; Saffaripour et al. 2014), which however, cannot be directly applied in aero engine simulations. Consequently, reduced mechanisms are needed, which usually have the drawback that they are valid for specific fuels or operating conditions, only (Xu et al. 2018; Yan et al. 2018). Fuel blends, modeled by a surrogate of representative species, are a further challenge (Stagni et al. 2014). Usually, only few species are considered in such surrogates, what results in a loss of information (Stagni et al. 2016; Zettervall et al. 2020). At the moment there is a gap between very detailed but computationally expensive models, and efficient ones, which are fuel or case dependent. This circumstance is addressed with the gas phase, $\hbox {PAH}$ and soot model presented in this work. Large parts of the $\hbox {PAH}$ and soot model are taken from previous works for ethylene (Blacha et al. 2011; Eberle et al. 2017). These models are further improved in this paper and an adaption to the new gas phase mechanism for complex fuels is achieved. Aim is to present a soot model with a broad range of applicability for different fuels and operating conditions, without adapting model parameters. The objective is not to outperform much more detailed approaches, but to enable accurate aero engine soot $\hbox {CFD}$ simulations within reasonable computational time.

The defined model is validated by test cases with various single component fuels and blends. A number of shock tube experiments and several atmospheric laminar sooting flames are investigated. Good agreement is achieved in most studied experiments. It should be mentioned that the validation of the $\hbox {PAH}$ model is achieved by indirect validation of gas phase chemistry and soot model. Either corresponding $\hbox {PAH}$ measurements are not provided for the experiments under consideration, or too few representatives of a $\hbox {PAH}$ section were measured, which may invalidate the comparison. The developed model is an extension of previous works to more complex fuels (Blacha et al. 2011; Eberle et al. 2017) with new features and improved reaction mechanisms.

2 Methodology

The $\hbox {German Aerospace Center (DLR)}$ combustion code ThetaCOM provides the framework for the implementation of the subsequently defined models. This unstructured finite-volume solver is optimized for reacting flows in the low Mach number regime. A $\hbox {FRC}$ approach is utilized to describe the chemical processes (Warnatz et al. 2006; Gerlinger 2005) which allows a correct treatment of gas phase, $\hbox {PAH}$ and soot interaction, and satisfies element and mass conservation. The $\hbox {FRC}$ approach covers inherently the relatively slow soot chemistry without additional modeling efforts.

Fig. 1

Discretization of the mass spectrum for gas phase, $\hbox {PAHs}$ and soot model. $\hbox {PAH}_{i}^{*}$ denotes a radical belonging to the respective $\hbox {PAH}$ section

Soot evolution is divided into three sub-models for gas phase, $\hbox {PAHs}$, and soot particles. A sketch of the corresponding ranges in the mass spectrum is given in Fig. 1. In this figure $\hbox {PAH}_i$ and $\hbox {PAH}_{i}^{*}$ are artificial species, which cover all stable and radical $\hbox {PAHs}$ in the corresponding mass range. Due to the logarithmic scaling (factor 2 for all sections) the model covers the range of soot mass present in most sooting flames with less than 30 soot sections (bins). In addition to acetylene, polyynes are known to have an influence on soot growth especially at higher temperatures (Krestinin 1994, 1998, 2000). Therefore, the models for gas phase, $\hbox {PAHs}$, and soot are supplemented by reactions for these additional species.

In this work the transport equations for soot particles are defined in analogy to gas phase species and $\hbox {PAHs}$ by an Eulerian approach. Molecular differential diffusion is considered for gas phase species and $\hbox {PAHs}$ to accurately predict the diffusive processes, which are highly relevant in laminar flames (Williams 2019; Kronenburg and Bilger 2001). Radiation is modeled assuming an optically thin medium, whereby the additional calculation effort is kept low, since reabsorption of emitted energy is neglected. This simplification is applicable to the laminar flames investigated in this work (Liu et al. 2002). In addition to soot, $\hbox {CO}_2$ and $\hbox {H}_2\hbox {O}$ are the most radiating species. The required coefficients are taken from the work of Mauss (1997).

2.1 Gas Phase Chemistry

The basis of the presented gas phase chemistry is a detailed chemical kinetic model for small hydrocarbons and $\hbox {PAHs}$ developed by Slavinskaya et al. (2012). From this work, the complete $\hbox {C}_{1}$-$\hbox {C}_{6}$ combustion is adopted, which also includes pathways for $\hbox {PAH}$ key intermediates, e.g., $\hbox {C}_{5}\hbox {H}_{5}$, $\hbox {C}_{5}\hbox {H}_{6}$, $\hbox {C}_6\hbox {H}_5\hbox {O}$, or $\hbox {C}_6\hbox {H}_5\hbox {OH}$. To allow the simulation of more complex fuels and even blends such as kerosene, the gas phase model is extended by additional reaction pathways for $\hbox {C}_{7}$-$\hbox {C}_{10}$ combustion, as described in Slavinskaya (2008). To account for polyynes in PAH and soot formation, additional species ($\hbox {C}_{4}\hbox {H}_{2}$, $\hbox {C}_{6}\hbox {H}_{2}$, $\hbox {C}_{8}\hbox {H}_{2}$, $\hbox {C}_{10}\hbox {H}_{2}$, and $\hbox {C}_{12}\hbox {H}_{2}$) and corresponding pathways are added according to Slavinskaya et al. (2012). For comparative purposes in Sect. 3.2, a reference mechanism is defined which additionally includes all stable and radical $\hbox {PAH}$ species from Slavinskaya et al. (2012).

The developed gas phase mechanism considers only species up to a molar mass of ${100}\,\hbox {g/mol}$. As a result, the mechanism does not include most of the $\hbox {PAHs}$, which are explicitly captured by the sectional approach (see Sect. 2.2). On the other hand, toluene ($\hbox {C}_{7}\hbox {H}_{8}$) which may represent an aromatic fuel component of kerosene, and benzene ($\hbox {A}_{1}$) are included in the gas phase mechanism. Fuel radicals, e.g., benzyl ($\hbox {C}_{7}\hbox {H}_{7}$) or phenyl ($\hbox {A}_1\hbox {m}$), are included in the kinetic scheme as well. To keep the computational effort low, some combustion intermediates (e.g., $\hbox {iC}_{3}\hbox {H}_{7}$, $\hbox {iC}_{4}\hbox {H}_{5}$ or $\hbox {iC}_{5}\hbox {H}_{9}$) are neglected. Some additional minor differences compared to the original mechanism of Slavinskaya et al. (2012) are a change in activation energy of reaction 457 (related to the kinetic scheme in Online Resource 2) which is increased from $22393\,\hbox {K}$ to $28000\,\hbox {K}$, while the reactions ${\hbox {C}_4\hbox {H}_4=2\hbox {C}_2\hbox {H}_2}$, ${\hbox {CH}_4 + \hbox {CH}_3=\hbox {C}_2\hbox {H}_5 + \hbox {H}_2}$ and ${\hbox {C}_4\hbox {H}_8=\hbox {C}_2\hbox {H}_4 + \hbox {C}_2\hbox {H}_4}$ have been removed. The final gas phase mechanism consists out of 112 species and 846 reactions, of which are 807 reversible and 39 irreversible. In agreement with the derived gas phase model, a blend of ${65}\,\hbox {vol-}\%$ n-decane, ${23}\,\hbox {vol-}\%$ iso-octane and ${12}\,\hbox {vol-}\%$ toluene is defined as a possible kerosene (Jet A-1) surrogate. The resulting mechanism is optimized in terms of numerical stability by replacing reactions 174, 201, 202, 236, 240, 258, 794, 798 and 825 (numbers correspond to the kinetic scheme in Online Resource 2) with substituents taken from the works of Sivaramakrishnan et al. (2006), Dagaut and Gaïl (2007), Blacha et al. (2011), Smith et al. (2012), Raj et al. (2012), Aghsaee et al. (2014). The highly instable reaction ${\hbox {CH}_3\hbox {OH} + \hbox {O}_2 = \hbox {CH}_2\hbox {OH} + \hbox {HO}_2}$ has been removed, and two reversible reaction equations (406 and 837) are replaced by irreversible ones. In this way a reduced and stability optimized kinetic gas phase mechanism is developed. The mechanism and corresponding thermo and transport data are provided in Online Resource 2. A detailed investigation by a series of ignition delay, pyrolysis and oxidation shock tube experiments follows in Sect. 3.2.

2.2 PAH Model

$\hbox {PAHs}$ are soot precursors on the way from gaseous species to the particle phase. Instead of considering all $\hbox {PAHs}$ included in the mechanism of Slavinskaya and Haidn (2011), most of these molecules are lumped into sections within a certain mass range. As a tradeoff between computing time and accuracy, three logarithmically scaled $\hbox {PAH}$ sections cover the range from ${100} \,\hbox {g/mol}$ to ${800}\,\hbox {g/mol}$ (Di Domenico et al. 2010) (see Fig. 1). Thus, even small $\hbox {PAHs}$, such as naphthalene, are captured by the first section. The upper bound of the last $\hbox {PAH}$ section coincides with the lower bound of the first soot section, which was chosen to obtain a similar mean molar mass and particle diameter to incipient soot particles (Abid et al. 2009). A constant particle number density is assumed for the intra-sectional $\hbox {PAH}$ mass distributions. The $\hbox {hydrogen/carbon (H/C)}$ ratio, as well as the thermodynamic and transport properties of each $\hbox {PAH}$ section are calculated from reference species by mass-weighted inter- and extrapolation (Yu et al. 2004; Richter et al. 2005).

The reversibility of the $\hbox {PAH}$ formation step is a key property addressed by this model. Therefore, each $\hbox {PAH}$ section is supplemented by a radical branch ($\hbox {PAH}^*$) to facilitate the implementation of the $\hbox {HACA}$ mechanism in a reversible formulation (Frenklach and Wang 1991). The properties of these $\hbox {PAH}$ radicals are taken from their stable counterparts. The elementary composition is reduced by one ${\hbox {H}}$ atom and the molar enthalpy is increased by a constant value of ${250}\,\hbox {kJ/mol}$ (Burcat and Ruscic 2005). Uncertainties in the thermodynamic properties of these pseudo species have no major impact since radical concentrations in typical gas mixtures are usually several orders of magnitude lower than for stable species.

A reaction scheme is defined to describe the interaction between gas phase and the first $\hbox {PAH}$ section. For this purpose, all reactions involving gas phase species and $\hbox {PAHs}$, in the range from ${100}\,\hbox {g/mol}$ to ${200}\,\hbox {g/mol}$, are taken from the reference mechanism. Among others, these include reactions with fuel radicals such as benzyl and phenyl. These reactions are separated by whether they contain a stable or a radical $\hbox {PAH}$ molecule, which is replaced by $\hbox {PAH}_{1}$ or $\hbox {PAH}^{*}_{1}$, respectively. Stoichiometric coefficients are modified to achieve element and mass conservation. The $\hbox {PAH}_{1}$/$\hbox {PAH}^{*}_{1}$ gas phase interaction reactions are summarized in Tables 2 and 3 in Online Resource 1.

In addition to the described gas phase interaction, the $\hbox {PAH}$ model also includes surface growth, collision, and oxidation. As already mentioned, surface growth is implemented according to the $\hbox {HACA}$ mechanism. First, $\hbox {PAH}$ radicals are formed either by unimolecular decomposition or by hydrogen abstraction with $\hbox {OH}$, ${\hbox {O}}$ or ${\hbox {H}}$. Evolving radicals then directly contribute to growth via carbon addition from acetylene and diacetylene. Here, diacetylene is considered as part of the polyyne extension. Surface growth rates of these species are obtained from reactions with naphthalene and styrene, given in the reference mechanism. This results in a system of reactions which is exemplarily shown for hydrogen abstraction via ${\hbox {H}}$ and carbon addition by $\hbox {C}_{2}\hbox {H}_{2}$, according to

$$\begin{aligned} \hbox {PAH}_i + \hbox {H}&\rightleftharpoons \nu ^{\prime \prime }_{{\mathrm{PAH}}_{i}^{*}}\ \hbox {PAH}_{i}^{*} + \nu _{{\mathrm{H}}_{2}}^{\prime \prime }\ \hbox {H}_{2} , \end{aligned}$$

(1)

$$\begin{aligned} \hbox {PAH}_{i}^{*} + \hbox {C}_{2}\hbox {H}_{2}&\rightarrow \nu _{{\mathrm{PAH}}_i}^{\prime \prime }\ \hbox {PAH}_i + \nu _{{\mathrm{PAH}_{i+1}}}^{\prime \prime }\ {\mathrm{PAH}_{i+1}} + \nu _{{\mathrm{H}}}^{\prime \prime }\ \hbox {H}\quad i \in [1,3] . \end{aligned}$$

(2)

Soot inception is neglected in Eq. (2) and thus, for $i=3$, the product $\nu _{{\mathrm{PAH}}_{4}}^{\prime \prime }\ \hbox {PAH}_{4}$, which is replaced by $\nu _{{\mathrm{SOOT}}_{1}}^{\prime \prime }\ \hbox {SOOT}_{1}$, is set to zero by the stoichiometric prefactor ($\nu _{{\mathrm{SOOT}}_{1}}^{\prime \prime }=0$). To ensure element and mass conservation $\nu _{{\mathrm{PAH}}_{3}}^{\prime \prime }$ and $\nu _{{\mathrm{H}}}^{\prime \prime }$ are modified. However, collisions involving the last $\hbox {PAH}$ or $\hbox {PAH}$ radical section contribute to soot inception. The transfer of mass, resulting from the collision of $\hbox {PAHs}$ with other molecules or $\hbox {PAH}$ sections, follows the $\hbox {simultaneous particle and molecule modeling (SPAMM)}$ approach described by Pope and Howard (1997). To additionally maintain the conservation of elements an iterative technique (Blacha et al. 2011) is used for the calculation of stoichiometric coefficients. If required, corresponding reactions are supplemented with $\hbox {H}$ or $\hbox {H}_{2}$. The reaction rates of these collision-driven processes are calculated according to the kinetic theory of gases by

$$\begin{aligned} k = 2.2\ N_A\ \gamma _{i,j}\ \beta _{i,j}\quad i,\ j \in [1,3] . \end{aligned}$$

(3)

Here, a van der Waals enhancement factor of 2.2 is applied (Harris and Kennedy 1988; Miller 1991) and $N_A$ represents the Avogadro constant. As only interactions with at least one $\hbox {PAH}$ radical are considered, a collision efficiency $\gamma _{i,j}$ of unity is chosen. The collision frequency $\beta _{i,j}$ corresponds to the number of molecular collisions in a volume per unit time (Kazakov and Frenklach 1998). The chosen calculation approach follows the definitions of Fuchs et al. (1965) and is valid over the entire Knudsen number range.

The introduced $\hbox {PAH}$ radicals allow the oxidation reactions and associated rate constants to be directly taken from the reference mechanism (Slavinskaya et al. 2012). ${\hbox {OH}}$, $\hbox {O}_{2}$ and $\hbox {O}$ are identified to be the most important $\hbox {PAH}$ oxidizers and are therefore considered. A summary of the above described $\hbox {PAH}$ formation and oxidation processes is given in Table 4 in Online Resource 1. Compared to Eberle et al. (2017), the $\hbox {PAH}$ model has been improved to satisfy element conservation. In addition, it uses slightly different surface growth rates which are taken from respective reactions in the mechanism of Slavinskaya et al. (2012).

2.3 Soot Model

The sectional soot model comprises 30 logarithmically scaled sections, which ensure a sufficiently fine resolution of the mass and particle size spectrum. The lower limit of the first section is set to ${800}\,\hbox {g/mol}$ and guarantees a seamless transition from $\hbox {PAHs}$ to soot. A mass growth factor of 2 is chosen. In this way incipient soot particles correspond with the first soot section (Abid et al. 2009), while the last section ($\hbox {SOOT}_{30}$) covers extremely large particles up to an average collision diameter of ${8.3}\,{\upmu \hbox {m}}$. A constant particle number density is assumed for the intra-sectional mass distribution, as in the $\hbox {PAH}$ model. Following Richter et al. (2005), a decreasing $\hbox {H/C}$ ratio with increasing soot mass is obtained, in the present case

$$\begin{aligned} H/C = 0.4405\ M_{{\mathrm{SOOT}_i}}^{-0.10524} , \end{aligned}$$

(4)

where $M_{{\mathrm{SOOT}_i}}$ ($\hbox {kg mol}^{-1}$) is the molar mass of the ith soot section (Blacha et al. 2011). Thermodynamic properties of soot are extremely rare in literature and are thus determined by a mass-weighted scaling from $\hbox {C}_{2}\hbox {H}_{2}$. The soot density $\rho _{SOOT}$ is according to Lindstedt (1994) set to ${1800}\,\hbox {kg/m}^3$. The primary particle diameter of spherical soot particles follows from

$$\begin{aligned} d_{p,{\mathrm{SOOT}_i}} = \root 3 \of {6\ M_{{\mathrm{SOOT}_i}} / (N_A\ \pi \ \rho _{SOOT})}. \end{aligned}$$

(5)

Thereby, the spherical diameter of the first soot section is approximately $1\,\hbox {nm}$. Based on the work of Eberle et al. (2017) and in agreement with experimental studies (Kazakov et al. 1995; Smooke et al. 2005; Saggese et al. 2015), a diameter of $14\,\hbox {nm}$ (applies to $\hbox {SOOT}_{11}$) is defined to initiate soot aggregation. It is assumed that soot particles up to this section have a spherical structure with increasing diameter. The diameter and the number of primary particles per soot aggregate ($n_{p,{\mathrm{SOOT}_i}}$) above this threshold value follow from

$$\begin{aligned} d_{p,{\mathrm{SOOT}_i}}&= \root 3 \of {2 / \mathcal{X}_{agg}}\ d_{p,{\mathrm{SOOT}_{i-1}}},\text { and} \end{aligned}$$

(6)

$$\begin{aligned} n_{p,{\mathrm{SOOT}_i}}&= \mathcal{X}_{agg}\ n_{p,{\mathrm{SOOT}_{i-1}}}. \end{aligned}$$

(7)

Here, the model parameter $\mathcal{X}_{agg}$ is introduced to achieve a blending between coalescence ($\mathcal{X}_{agg}=1$) and aggregation ($\mathcal{X}_{agg}=2$). Based on the simulations in Sect. 3.4, best results are achieved with a value of 1.5 which is a compromise between pure coalescence and aggregation.

Similar to gas phase species and $\hbox {PAHs}$, each soot section is treated by a transport equation. Due to high aerosol Schmidt numbers, molecular diffusion can be neglected for soot particles. However, thermophoresis causes a drift of soot particles which scales with the temperature gradient and has to be considered. Using Einstein’s summation convention, the transport equation for the mass fraction of a soot section is given by

$$\begin{aligned} \frac{\partial }{\partial t} (\rho Y_{{\mathrm{SOOT}_i}}) + \frac{\partial }{\partial x_j} (\rho u_j Y_{{\mathrm{SOOT}_i}})&= \frac{\partial }{\partial x_j} \left( C_{th} \frac{\mu Y_{{\mathrm{SOOT}_i}}}{T}\frac{\partial T}{\partial x_j}\right) + \omega _{Y_{{\mathrm{SOOT}_i}}}. \end{aligned}$$

(8)

Therein $\rho$ denotes the density, $Y_{{\mathrm{SOOT}_i}}$ the mass fraction of the ith soot section, $u_j$ the velocity, $C_{th}$ the thermophoretic coefficient, $\mu$ the viscosity, T the temperature and $\omega _{Y_{{\mathrm{SOOT}_i}}}$ the chemical source term. Under the condition of a free molecular flow, the thermophoretic coefficient for soot is set to 0.55 (Messerer et al. 2003).

The soot chemistry includes the processes of surface growth, $\hbox {PAH}$-soot collision, agglomeration, and oxidation. Soot inception is already implemented in the $\hbox {PAH}$ collision model. Soot surface growth reactions are modeled in a simplified way based on the $\hbox {HACA}$ mechanism (Frenklach and Wang 1991). Because the number of soot sections is much higher than for the $\hbox {PAHs}$, soot radicals are neglected to keep the computational effort low. Instead, a quasi-steady state is assumed for these intermediate species and the reaction rate for acetylene or polyyne addition is given by

$$\begin{aligned} k_{{\mathrm{SOOT}_i},growth} = k_{C_xH_y}\ \mathcal{X} \ r_{\upchi }\ A_{{\mathrm{SOOT}_i}}\ \alpha . \end{aligned}$$

(9)

Here, $k_{C_xH_y}$ is the rate constant for the primary soot growth contributor $\hbox {C}_{2}\hbox {H}_{2}$ and additional polyynes, taken from Naydenova et al. (2004) and Wen et al. (2006). Following Frenklach and Wang (1991), $\mathcal{X}$ denotes the density of stable $\hbox {H/C}$ sites and is set to ${2.32\times 10^{-19}}\,{1/\hbox {m}^2}$. The factor $r_{\upchi }$ describes the ratio of reactive to stable sites and requires calibration. A value of ${1.7\times 10^{-3}}$ was found to be best suited for all investigated cases and is well within the range of literature values (Woods and Haynes 1991; Eberle et al. 2017). It should be noted, that the parameters $\mathcal{X}$ and $r_{\upchi }$ are the same for all soot sections. The soot surface per unit volume for primary particles is defined by

$$\begin{aligned} A_{{\mathrm{SOOT}_i}} = \pi \ d_{p,{\mathrm{SOOT}_i}}^2\ n_{p,{\mathrm{SOOT}_i}}\ C_{{\mathrm{SOOT}_i}} , \end{aligned}$$

(10)

with the soot concentration $C_{{\mathrm{SOOT}_i}}$ ($\hbox {mol m}^{-3}$). Finally, the parameter $\alpha$ depicts a further temperature dependency encountered in these reactions. To preserve the physical range of $0< \alpha < 1$, and to maintain the Arrhenius like form of the chemical reactions, $\alpha$ is defined as

$$\begin{aligned} \alpha (T,c_{\alpha },T_{\alpha },n_{\alpha }) = c_{\alpha }\ (T_{\alpha } / T)^{n_{\alpha }}\ \exp {(n_{\alpha }\ [1 - T_{\alpha } / T])} . \end{aligned}$$

(11)

The coefficients $c_{\alpha }=1$, $T_{\alpha }=1800$ and $n_{\alpha }=40$ are taken from the work of Eberle (2020). As for the $\hbox {PAHs}$ in Eq. (3), reaction rates for soot collisions are obtained from the kinetic theory of gases. Stoichiometric coefficients are determined in analogy to Sect. 2.2 with the $\hbox {SPAMM}$ approach. Based on the investigations of Blacha et al. (2011), the collision efficiency is set to 0.3 for interactions with $\hbox {PAHs}$ and $\hbox {PAH}$ radicals. An efficiency of unity is assumed for soot agglomeration reactions. Since soot aggregates differ in their properties to spherical soot particles, the collision diameter $d_{c,{\mathrm{SOOT}_i}}$ is calculated in terms of a gyration radius $R_{gyr,{\mathrm{SOOT}_i}}$ as proposed by Koeylue et al. (1995)

$$\begin{aligned} d_{c,{\mathrm{SOOT}_i}}&= 2\ R_{gyr,{\mathrm{SOOT}_i}}\ \sqrt{(D_f + 2) / D_f} , \end{aligned}$$

(12)

$$\begin{aligned} R_{gyr,{\mathrm{SOOT}_i}}&= (n_{p,{\mathrm{SOOT}_i}} / k_g)^{1 / D_f}\ d_{p,{\mathrm{SOOT}_i}} / 2 . \end{aligned}$$

(13)

In general, the fractal dimension $D_f$ is 3 for spherical, 2 for planar, and 1 for linear structures. A relatively universal value of $D_f={1.7\pm 0.15}$ for soot aggregates with $24\,\hbox {nm}< d_{p,{\mathrm{SOOT}_i}} < 52\,\hbox {nm}$ has been found for a large number of fuels, investigated in various laminar and turbulent flames (Koeylue et al. 1995). A molecule structure dependent transition of $D_f$ is adopted from Rosner and Pyykönen (2002), which originally characterizes the average fractal dimension of a homogeneous mixture. $D_f$ is obtained considering a non-sphericity parameter $\xi = \root 3 \of {n_{p,{\mathrm{SOOT}_i}}} - 1$ as follows

$$\begin{aligned} D_f = {\left\{ \begin{array}{ll} {1.7} &{} \text {if } \xi \ge 1 ,\\ {1.7} + {1.3} \cos {(0.5\ \pi \ \xi )} &{} \text {otherwise} . \end{array}\right. } \end{aligned}$$

(14)

The dependency on $\xi$ ensures a smooth transition from primary particles ($D_f = 3$) to soot aggregates ($D_f = 1.7$). Furthermore, the required lacunarity is calculated according to Rosner and Pyykönen (2002) by

$$\begin{aligned} k_g = \sqrt{3} + 0.21 (D_f - 1) + 2.29 / 4 (D_f - 1)(3 - D_f) . \end{aligned}$$

(15)

Soot oxidation is dominated by the species $\hbox {O}_{2}$ and $\hbox {OH}$, with rates taken from Fenimore and Jones (1967), and Appel et al. (2000), respectively. While the oxidation with $\hbox {O}_{2}$ is implemented following Eq. (9), again with $r_{\upchi }={1.7\times 10^{-3}}$ and identical $\alpha$ as for soot surface growth, the oxidation via $\hbox {OH}$ follows from

$$\begin{aligned} k_{{\mathrm{SOOT}_i},ox} = \eta \ k_{OH}\ N_A\ A_{{\mathrm{SOOT}_i}}\ \alpha . \end{aligned}$$

(16)

Based on established literature values (Neoh et al. 1981; Edwards et al. 2014), the efficiency $\eta$ is set to 0.13. Moreover, Liu et al. (2003) observed a decrease in $\hbox {OH}$ oxidation efficiency below a threshold temperature. To implement such a behavior while maintaining an Arrhenius like form, a combination of three $\hbox {OH}$ oxidation reactions is introduced. Each reaction is multiplied with a Gaussian pulse shaped function $\alpha$, as defined in Eq. (11), to achieve the observed decrease in oxidation efficiency by superposition. The coefficients for these additional functions are obtained by a least-square optimization to $\alpha _1 = \alpha (T,0.88,1826.58,70)$, $\alpha _2 = \alpha (T,0.58,2293.41,70)$ and $\alpha _3 = \alpha (T,0.94,3055.95,28.27)$. The formation and oxidation processes are summarized in Table 5 in Online Resource 1.

3 Results and Discussion

The presented gas phase, $\hbox {PAH}$ and soot model is validated using a large number of experiments with different operating conditions and fuels. Only a limited subset is presented here, which aims at the application in aero engines. Thereby, model features such as the polyyne extension are examined, and the choice of model parameters is discussed. Apart from the following investigation of the gas phase mechanism and the soot model, an individual validation of the $\hbox {PAHs}$ is not feasible. Either measurements for these species are not provided at all in the experiments, or data for few $\hbox {PAH}$ molecules are measured only, which is not enough, to enable a meaningful comparison with the lumped $\hbox {PAH}$ sections. Therefore, the $\hbox {PAH}$ model is validated here indirectly only, using soot data which is linked by the $\hbox {PAHs}$ to the gas phase.

First, shock tube experiments are used to examine ignition delay times, species profiles and soot yields. Fuel burn-out and soot oxidation can be investigated by these experiments to a limited extent only. However, the oxidation of soot in $\hbox {Rich-Burn, Quick-Mix, Lean-Burn (RQL)}$ combustion chambers is promoted by the addition of secondary air downstream of the main reaction zone. A similar behavior is observed in gaseous coflow diffusion flames, where ambient air surrounds the combustion zone. These flames are therefore well suited for a further validation and offer the advantage that no uncertainties are introduced by the modeling of liquid phase and spray. To take also premixing into account, a partially premixed ethylene flame (McEnally and Pfefferle 2000) is simulated and discussed in Online Resource 1.

3.1 Ignition Delay Times

Fig. 2

Measured and calculated ignition delay times: a ethylene and benzene, and b toluene, iso-octane and n-decane. Open symbols: experimental data (Horning 2001; Burcat et al. 1985); lines with filled symbols: calculated data

The experimental characterization of ignition delay times allows a first validation of the gas kinetic scheme. The compared ignition delay time is defined here to be the time difference between the point of initialization and maximum temperature gradient. Results for ethylene ($\Phi =1$, $p=1\,\hbox {bar}$), benzene ($\Phi =2$, $p={2}\,\hbox {bar}$ to ${3}\,\hbox {bar}$), toluene ($\Phi =1$, $p={2}\,\hbox {bar}$ to ${3}\,\hbox {bar}$), iso-octane ($\Phi =2$, $p=1\,\hbox {bar}$) and n-decane ($\Phi =1$, $p=1\,\hbox {bar}$) are shown in Fig. 2a, b together with experimental data taken from Horning (2001) and Burcat et al. (1985). As apparent, the gas phase mechanism achieves an excellent agreement with the experiments for ethylene, toluene and n-decane. It is observed that the ignition delay time for benzene is slightly underestimated at higher temperatures and for iso-octane at lower temperatures. Iso-octane is especially important because it is part of the chosen kerosene surrogate. In contrast to n-decane, however, the amount of iso-octane in the kerosene surrogate in this work is minor, and therefore a significant impact on the ignition delay time is not expected.

3.2 Shock Tube Species Profiles

Fig. 3

Measured and calculated species profiles for the oxidation of kerosene at three stoichiometric ratios: a acetylene ($\hbox {C}_{2}\hbox {H}_{2}$), b ethylene ($\hbox {C}_{2}\hbox {H}_{4}$), and c benzene ($\hbox {C}_{6}\hbox {H}_{6}$). Open symbols: experimental data (Malewicki et al. 2013); solid lines with filled symbols: calculated data; dash-dotted lines with filled symbols: calculated data from the reference gas phase mechanism (ref)

Species profile measurements behind reflected shockwaves are examined to further validate the gas phase mechanism. Experiments for methane, ethylene, acetylene, benzene, n-heptane, and components of the kerosene surrogate are simulated under various operating conditions. Here, an experiment for the oxidation of kerosene is chosen which addresses the purpose of this work, an application of the gas phase mechanism in aero engines. A good prediction of combustion products as well as of soot intermediates is essential.

Malewicki et al. (2013) investigated three kerosene–air equivalence ratios ($\Phi =0.46, 0.98\, {\text{and}}\, 1.86$) at pressures ranging from ${16.3}\,\hbox {bar}$ to ${28.7}\,\hbox {bar}$ and total reaction times from ${1.24}\,\hbox {ms}$ to ${3.32}\,\hbox {ms}$. The important $\hbox {PAH}$ and soot surface growth species $\hbox {C}_{2}\hbox {H}_{2}$, and the fuel decomposition product $\hbox {C}_{2}\hbox {H}_{4}$ are in excellent agreement with the measured values as shown in Fig. 3a, b. A slightly growing deviation from the experiment with increasing $\Phi$ is observed for benzene in Fig. 3c. However, the tendency with increasing $\Phi$ is well reproduced. The achieved soot yields are similar to those of Malewicki (2012). It also has to be mentioned that there is an uncertainty due to the composition of experimentally investigated and simulated kerosene fuel.

Figure 3 additionally shows for $\Phi =0.46$ a comparison of the newly developed gas phase, $\hbox {PAH}$ and soot model with the reference mechanism (Slavinskaya and Haidn 2011), from which large parts of this model are derived. The reference mechanism includes 243 species and 1634 reactions, treats $\hbox {PAHs}$ as physical species and neglects soot. The derived mechanism provides a considerably better agreement with the experiment (compare with dash-dotted lines). Since both models have a very similar gas phase mechanism, the extended mass range of the $\hbox {PAH}$ model probably causes the strong decrease of $\hbox {C}_{2}\hbox {H}_{2}$, which in turn affects the formation of $\hbox {C}_{2}\hbox {H}_{4}$ and benzene. In particular, the improved prediction of benzene, which dominates the interaction between the first $\hbox {PAH}$ section and gas phase in this model, confirms a high accuracy of the defined mechanism for the presented operating conditions. Corresponding simulations were also carried out without activated soot model to exclude influences on the results. Due to the very low soot concentrations in this case the model achieves identical results to the model which considers soot as shown in Fig. 3. Thereby, the reference mechanism leads to a 1.5 times higher computational effort than for the model with soot. Without soot model, the factor is 2.6.

3.3 Shock Tube Soot Yields

Fig. 4

Measured and calculated soot yields in the pyrolysis and oxidation of various hydrocarbon fuels: a methane, b ethylene, c, d benzene, e toluene, f toluene/iso-octane, g, h n-heptane, and i n-decane. Open symbols: experimental data (Kellerer et al. 1996; Vlasov and Warnatz 2002; Agafonov et al. 2011, 2007; Alexiou and Williams 1995; Hong et al. 2007; Bockhorn 2013); solid lines with filled symbols: calculated data; dotted lines with filled symbols: calculated data without polyynes (w/o polyynes)

Table 1 Investigated fuels, boundary conditions and literature references for the chosen shock tube soot yield experiments

The soot emissions of various fuels at different operating conditions are also investigated for shock tube experiments. Here, the soot yield represents the concentration of carbon bounded in soot (in the simulation of all soot sections) compared to the total carbon concentration. Methane and ethylene are considered together with more complex fuels to demonstrate a wide scope of application of the presented soot model. A summary of investigated fuels with applied boundary conditions is given in Table 1. In addition, the influence of the polyyne extension is demonstrated in experiments covering the high temperature range (dotted lines).

The oxidation of methane is investigated at an equivalence ratio of $\Phi =5$ (Kellerer et al. 1996). As can be seen in Fig. 4a, the position and value of the calculated maximum soot yield are in good agreement with the experiment. In the section of the rise in soot yield, soot inception and growth are dominant. With increasing temperature, the thermal instability of the $\hbox {PAH}$ molecules become predominant, what causes the following decrease of soot yield. This behavior is captured in the simulation by the reversible formulation of the $\hbox {PAH}$ model, which reproduces the measured distribution very well. The following pyrolysis experiment with ethylene (see Fig. 4b) uses two different initial carbon concentrations (Vlasov and Warnatz 2002). Good agreement is obtained for the experiment with the lower carbon content. For the higher initial value, no measured data is available in the range from $1800\,\hbox {K}$ to $2100\,\hbox {K}$, and therefore no comparison to the experiment is possible in this relevant range of temperature. Despite this gap in the measured values, the model reflects the experimentally encountered trend with increasing carbon content very well.

In the next experiment (Vlasov and Warnatz 2002), benzene is investigated at identical operating conditions as before (illustrated in Fig. 4c). Position and value of the maximum soot yield, as well as the distribution and trends, agree well with the measured values for both cases. Compared to the simulations in Vlasov and Warnatz (2002) with a mechanism of similar size, the decrease of soot yield at higher temperatures is predicted considerably better, while deviations from the experiment at lower temperatures are comparable. In addition, the influence of polyynes on soot growth is highlighted (compare to dotted lines for model without polyynes). Up to a temperature of $2000\,\hbox {K}$, the model with polyynes agrees better with the experiment. However, as the temperature increases further, the measured soot yields are overestimated if polyynes are considered. A similar experiment (see Fig. 4d) investigates benzene for a reduced pressure (Agafonov et al. 2011). Compared to the experiments, the formation of soot is initiated at slightly lower temperatures. The model provides a good agreement with the experiment at temperatures above $1800\,\hbox {K}$ which is due to the consideration of polyynes. In the simulation with the highest $\hbox {C}_{6}\hbox {H}_{6}$ concentration but without polyynes, the soot yield in the high temperature region is too low.

In the oxidation experiment of toluene (Agafonov et al. 2007) (see Fig. 4e), the soot yield is overestimated by a factor of 1.6 in the simulation. Nevertheless, the model predicts the trend with increasing oxygen content correctly and reproduces the distributions with varying temperature very well. For this test case similar results are obtained as in the simulation of Agafonov et al. (2007) (210 species and 2250 reactions), however with a considerably smaller computational effort. As before the soot yield is overpredicted if polyynes are included in the high temperature range, especially in the rich case. The pyrolysis of iso-octane is examined in a mixture with toluene (Alexiou and Williams 1995) (see Fig. 4f). The maximum value of the soot yield as well as the position are properly predicted by the model. In comparison to the experiment the distribution of the soot yield extends over a wider range. This deviation is probably caused by the observed overestimation of soot in the toluene experiment (compare Fig. 4e) but is still in a range to be accepted.

Two further experiments provide insight into the oxidation of n-heptane (Agafonov et al. 2007; Hong et al. 2007). Due to its relevance to the automotive industry more experimental data is available for this fuel than for n-decane. However, given the similarity between both fuels, these experiments provide much information on the sooting tendency of a major kerosene component. The first experiment analyses the pressure dependency of soot evolution for $\Phi =5$. The second one examines soot trends at different reaction times for the same equivalence ratio. The results are visualized in Fig. 4g, h. In both experiments, the model can predict the positions and values of maximum soot yield, as well as distributions and trends accurately. Compared to soot yields predicted in Agafonov et al. (2007) with a more detailed gas phase chemistry, a considerably better agreement to measured values is achieved. A single experiment for n-decane has been found in literature only (see Fig. 4i). The conditions of this oxidation experiment are $\Phi =5$, a pressure of $0.5\,\hbox {bar}$ and a reaction time of ${1}\,{\hbox {ms}}$ (Bockhorn 2013). Except for the values, the distribution and the position of the soot yield maximum is reproduced precisely. However, the maximum value is by a factor of 4 higher than in the experiment.

Some of the cited experimental studies (Vlasov and Warnatz 2002; Agafonov et al. 2011, 2007) have been supplemented by numerical simulations, too. In these studies, either very detailed kinetic models were used and/or adjusted reaction rates have been proposed to achieve best results. Even in comparison with these simulations, the presented soot model provides good results, with a relatively small gas phase chemistry, and reproduces the observed trends very well with a single set of model constants. With respect to the polyynes the findings are inconclusive and due to existing measurement uncertainties, no distinct statement can be made here. However, it is demonstrated, that polyynes can be relevant at high temperatures above $2000\,\hbox {K}$, a range reached in aero engines. In this temperature range the soot yield is always increased if polyynes are considered.

3.4 Laminar Ethylene Diffusion Flame

Table 2 Boundary conditions for the laminar non-smoking and smoking ethylene diffusion flame (Santoro et al. 1987)

In a next step, the fuel flexibility as well as the oxidation properties of the $\hbox {PAH}$ and soot model are investigated. Test case is an atmospheric laminar sooting ethylene diffusion flame of Santoro et al. (1987). The experiment is operated with varying flow velocities of fuel and ambient air, to achieve either a non-smoking (F2) or a smoking (F4) flame. The term non-smoking refers to a soot emission below the visual and experimental detection limit at the combustion chamber exit. This test case has been studied by many groups, but most of them examined either the smoking or the non-smoking flame. Thus, it is not shown whether the transition from a non-smoking to a smoking flame is reproduced by their models (D’Anna et al. 2007; Zhang et al. 2009; Dworkin et al. 2011). In this work, corresponding simulations are performed on an axisymmetric grid (${5}^{\circ }$ wedge) with a height of $160\,\hbox {mm}$ and approximately 30000 volumes. The radial dimensions as well as the chosen boundary conditions are summarized in Table 2. With an offset of $4\,\hbox {mm}$ to the combustion chamber, a fully developed laminar inflow profile is predefined at the inlet of the fuel feed tube.

Fig. 5

Calculated data for a non-smoking (left) and smoking (right) laminar ethylene diffusion flame: a temperature, b soot volume fraction, and c $\hbox {OH}$ mass fraction. Superimposed are iso-lines of: significant heat radiation source term $S_r$ (solid), and mass fraction $Y_{{\mathrm{C}_2\mathrm{H}_x}}={1}$‰ of ethylene $Y_{{\mathrm{C}}_{2}{\mathrm{H}}_{4}}$ (dashed) and acetylene $Y_{{\mathrm{C}}_{2}{\mathrm{H}}_{2}}$ (dotted) in (a); soot source terms $\omega _S={0.1}\,{\hbox {g}/(\hbox {m}^3\hbox {s})}$ of inception $\omega _S^i$ (dotted), $\hbox {PAH}$-soot collision $\omega _S^p$ (dashed) and $\hbox {C}_{2}\hbox {H}_{2}$/polyyne addition $\omega _S^g$ (solid) in (b); two temperatures of $1600\,\hbox {K}$ (dashed) and $1800\,\hbox {K}$ (solid), and $\hbox {OH}$ oxidation soot source term $\omega _S^o={-0.1}\,{\hbox {g}/(\hbox {m}^3\hbox {s})}$ (solid, white) in (c)

Figure 5 shows predicted temperature, soot volume fraction and $\hbox {OH}$ mass fraction contour plots. In case of the smoking flame (always on the right), areas of high fuel mass fraction (see dashed line for $Y_{{\mathrm{C}}_{2}{\mathrm{H}}_{4}}$) extend further downstream. This is caused by a higher inflow velocity and is typical for laminar flames. Thus, in the lower part up to $h\approx 70\,\hbox {mm}$, the smoking flame is stretched in axial direction compared to the non-smoking one (see also temperature iso-lines in Fig. 5c). Due to high soot concentrations, however, temperature distributions get more complex further downstream. Soot radiation (see solid line in Fig. 5a) causes significant energy losses in both flames. While the region with high soot concentration and subsequently high radiative losses extends over the complete height in case of the smoking flame, it is limited to $h<70\,\hbox {mm}$ for the non-smoking one. Consequently, heat production exceeds heat losses further downstream, and a second high temperature region is obtained at $h\approx 70\,\hbox {mm}$ (see Fig. 5c) in the non-smoking case.

Figure 5b shows the predicted soot volume fraction distributions, superimposed with iso-lines of notable soot production terms for inception $\omega _S^i$ (dotted lines), collision with $\hbox {PAHs}$ $\omega _S^p$ (dashed lines) and $\hbox {C}_{2}\hbox {H}_{2}$/polyyne addition $\omega _S^g$ (solid lines). For both flames, the region of highest soot inception (outer dotted region in Fig. 5b) is located directly downstream of the injector tube wall, well within the inner, fuel rich side of the flame (see also dotted line for acetylene mass fraction $Y_{{\mathrm{C}}_{2}{\mathrm{H}}_{2}}={1}$‰ in Fig. 5a) and on the axis of the burner. The areas of significant source terms due to $\hbox {PAH}$-soot collision and $\hbox {C}_{2}\hbox {H}_{2}$ and polyyne addition are nearly identical. Minor differences appear close to the burner exit only, where $\omega _S^p$ increases a little bit faster. A comparison between the smoking and the non-smoking flame shows, that in the first case, due to the higher inflow velocity, soot formation starts with a slight offset but extends significantly further downstream. Especially, in the second part of the smoking flame, soot formation is not completely compensated by oxidation (region of $\hbox {OH}$ oxidation, see Fig. 5c).

Fig. 6

Measured and calculated profiles in an atmospheric laminar ethylene diffusion flame: a radial temperature (non-smoking), b radial $\hbox {OH}$ mole fractions (non-smoking), and c radially integrated soot volume fraction (smoking and non-smoking). Open symbols: experimental data (Santoro et al. 1987); solid lines with filled symbols: calculated data

Fig. 7

Measured and calculated radial profiles in an atmospheric laminar ethylene diffusion flame: a axial velocity (non-smoking), b axial velocity (smoking), and c ethylene mole fraction (non-smoking). Open symbols: experimental data (Santoro et al. 1987, 1983); solid lines with filled symbols: calculated data

The experimental study provides measured temperature data for the non-smoking flame, as shown in Fig. 6a. A deviation between simulation and experiment is observed near the burner inlet ($h=3\,\hbox {mm}$) on the axis. As there is no significant heat radiation in this area (see Fig. 5a), the largest uncertainty here is the thermal inflow boundary condition. The simulation uses the experimentally specified inflow boundary temperature of $300\,\hbox {K}$. $3\,\hbox {mm}$ further downstream $530\,\hbox {K}$ are obtained in the simulation mainly due to released energy from chemical conversion and heat conduction. With an excellent agreement between the calculated and measured velocity profiles (compare Fig. 7a, b) it remains unclear, why $750\,\hbox {K}$ are measured in the experiment after such a short distance. A possible cause could be the heating of the fuel tube in the experiment and subsequently a higher fuel inflow temperature, which is not considered in the simulation. The calculated maximum soot volume fractions at a height of $50\,\hbox {mm}$, are a factor of 1.3 above the measured values, causing higher soot radiation and lower temperatures close to the axis ($\text {radius} < 5\,\hbox {mm}$). Apart from this deviation, however, predicted temperatures are in a good agreement with the experiment. Figure 6b shows calculated and measured $\hbox {OH}$ mole fractions, an important quantity for soot oxidation. The radial trends and values are very well reproduced at both heights. While $\hbox {OH}$, $\hbox {O}_{2}$ and $\hbox {O}$ are species which are involved in soot oxidation, this experiment is dominated by $\hbox {OH}$ oxidation. The radially integrated soot volume fractions along the burner height are shown in Fig. 6c for the non-smoking and smoking flame. Compared to the experimental measurements, soot formation is initiated by the model too early. This leads to an overestimation of the integrated soot volume fraction by a factor of up to 2. At a height of $50\,\hbox {mm}$ the maximum values are very well predicted. As observed in the experiment, soot oxidation in the model also initiates here. In both cases, the model is able to reflect the trend over the burner height. In the non-smoking case predictions above a height of $85\,\hbox {mm}$ are below the last point measured and thus, also confirm the previously achieved visual classification between non-smoking and smoking (Fig. 5b).

Fig. 8

Measured and calculated radial profiles in an atmospheric laminar ethylene diffusion flame: a acetylene mole fraction (non-smoking), b soot volume fraction (non-smoking), and c soot volume fraction (smoking). Open symbols: experimental data (Santoro et al. 1983); solid lines with filled symbols: calculated data

Despite the underestimated temperatures for the non-smoking flame near the burner exit the consumption of ethylene is slightly higher according to the gas phase mechanism as in the experiment (see Fig. 7c). This in turn leads to an increased formation of acetylene at a height of $7\,\hbox {mm}$ as can be seen in Fig. 8a. The resulting premature formation of soot leads to the overestimation of the radial soot volume fraction by a factor of 2.5 at a height of $15\,\hbox {mm}$ in Fig. 8b. Moreover, a radial offset of $0.7\,\hbox {mm}$ between measured and predicted profile peak can be observed. The low temperatures on the axis and the overestimated formation of soot may lead to the acetylene profile at $20\,\hbox {mm}$ in Fig. 8a. At a height of $50\,\hbox {mm}$, the radial shift and the peak soot volume fraction deviation factor is $0.6\,\hbox {mm}$ and 1.3, respectively. The model yields nearly identical temperature and species profiles compared to other numerical studies (Dworkin et al. 2011; D’Anna et al. 2007). The radial peak offset, as well as the underestimation of soot volume fraction on the burner axis is also predicted in many other works (D’Anna et al. 2007; Dworkin et al. 2011; Zhang et al. 2009; Liu et al. 2003). The overestimation of the radial peak is found in D’Anna et al. (2007) as well. Predicted radial profiles for the smoking flame are compared with the measurements in Fig. 8c. The radial trend as well as the position of the maximum peak agree very well with the experiment at a height of $40\,\hbox {mm}$. However, at $70\,\hbox {mm}$, a similar deviation arises as found for the non-smoking flame. The maximum soot peak is slightly shifted towards the fuel-rich side and the soot volume fraction measured on the axis is underestimated by a factor of 3 (${1.5}\,\hbox {ppm}$). In comparison to other numerical investigations the model achieves good results for the smoking flame (Liu et al. 2003; Khosousi and Dworkin 2015; Akridis and Rigopoulos 2016). In fact, no model has been found that correctly predicts the radial trends, peak and axis values for the soot volume fraction at the same time.

3.5 Laminar n-Decane and Kerosene Diffusion Flame

Table 3 Boundary conditions for the laminar n-decane and kerosene diffusion flame (Saffaripour et al. 2014, 2011)

In the works of Saffaripour et al. (2011, 2013, 2014) the sooting tendencies of laminar atmospheric n-decane and kerosene diffusion flames are investigated experimentally and numerically. The same burner is used in both experiments which are operated at different inflow velocities for fuel and ambient air, as specified in Table 3. The simulations are performed on an axisymmetric grid (${5}^{\circ }$ wedge) with a height of $120\,\hbox {mm}$ and approximately 25000 volumes. Sufficient spatial resolution for the chosen grid has been verified using a grid study comprising four levels of refinement. In comparison to the simulations with a higher resolution, this grid achieves almost identical results with a considerably less computational effort. As for the ethylene flame in Sect. 3.4, a fully development laminar inflow profile is specified for the fuel feed.

Fig. 9

Calculated data for an atmospheric laminar n-decane (right) and kerosene (left) diffusion flame: a temperature, b soot volume fraction, and c $\hbox {C}_{2}\hbox {H}_{2}$ mass fraction. Superimposed are iso-lines of: fuel mass fractions Y_f,1‰ (solid, white) and $Y_{f,{1}\,{\mathrm{ppm}}}$ (dashed, white), and mass fraction ratio $Z_i=Y_i/Y_{i,\text {inflow}}=1\%$ of iso-octane $Z_{{\mathrm{iC}}_{8}{\mathrm{H}}_{18}}$ (solid), n-decane $Z_{{\mathrm{nC}}_{10}{\mathrm{H}}_{22}}$ (dashed) and toluene $Z_{{\mathrm{C}}_{7}{\mathrm{H}}_{8}}$ (dotted) in (a); mass fractions $Y={1}\,\hbox {ppm}$ of incipient soot $Y_{{\mathrm{SOOT}_{1}}}$ (solid), first soot aggregates $Y_{\mathrm{SOOT}_{11}}$ (dotted) and soot particles of largest mass $Y_{\mathrm{SOOT}_{30}}$ (dashed) in (b); soot source terms $\omega _S={1}\,{\hbox {g}/(\hbox {m}^3\hbox {s})}$ of inception $\omega _S^i$ (dotted), $\hbox {C}_{2}\hbox {H}_{2}$ and polyyne addition $\omega _s^g$ (solid), $\hbox {PAH}$-soot collision $\omega _S^p$ (dashed), and $\hbox {OH}$ oxidation $\omega _S^o={-1}\,{\hbox {g}/(\hbox {m}^3\hbox {s})}$ (dash-dotted) in (c)

The contour plots of calculated temperature, soot volume fraction and $\hbox {C}_{2}\hbox {H}_{2}$ mass fraction are visualized in Fig. 9. The plots for n-decane are always shown on the left, and for kerosene on the right side. In case of n-decane (see Fig. 9a), iso-lines of constant fuel mass fraction $Y_{f,{1}}$_‰ and temperature reach further downstream than for the kerosene flame. This is expected due to a higher inflow velocity but could also be affected by chemical kinetics. Moreover, it can be observed that the fuel is consumed much faster in the case of the n-decane flame (compare Y_f,1‰ and $Y_{f,{1}\,{\mathrm{ppm}}}$ between both flames). The decay of the three kerosene surrogate components is shown in Fig. 9a on the right side. In contrast to n-decane $Z_{\mathrm{nC}_{10}\mathrm{H}_{22}}$ (dashed line), iso-octane $Z_{\mathrm{iC}_{8}\mathrm{H}_{18}}$ (solid line) is consumed faster, and toluene $Z_{\mathrm{C}_{7}\mathrm{H}_{8}}$ (dotted line) persists still longer and in areas of higher temperature. $Z_i$ is the ratio of the local species mass fraction related to its inflow value. In the n-decane flame less soot is obtained than in the kerosene flame, but both show a similar non-smoking behavior. However, the comparison of individual soot sections in Fig. 9b, shows that a small amount of largest soot particles $Y_{\mathrm{SOOT}_{30}}$ (dashed lines) is not completely oxidized in the kerosene flame. Iso-lines of incipient soot $Y_{{\mathrm{SOOT}_{1}}}$ (solid lines) and first soot aggregates $Y_{\mathrm{SOOT}_{11}}$ (dotted lines) are very similar in both flames, except for their stretching in flow direction. For the n-decane flame, there are almost no overlapping areas between these three shown soot sections. Figure 9c shows the predicted $\hbox {C}_{2}\hbox {H}_{2}$ mass fraction. Here, the superimposed soot source terms of inception $\omega _S^i$ (dotted lines), $\hbox {C}_{2}\hbox {H}_{2}$ and polyyne addition $\omega _s^g$ (solid lines) and $\hbox {PAH}$-soot collision $\omega _S^p$ (dashed lines) are found as expected in the presence of $\hbox {C}_{2}\hbox {H}_{2}$. Again, oxidation by $\hbox {OH}$ $\omega _S^o$ (dash-dotted lines) is dominant and soot oxidation by $\hbox {O}_{2}$ is minor in this case.

Fig. 10

Measured and calculated profiles in an atmospheric laminar n-decane diffusion flame: a axial soot volume fraction, b radial soot volume fraction, and c radial temperature. Open symbols: experimental data (Saffaripour et al. 2014, 2011); solid lines with filled symbols: calculated data

Axial and radial soot volume fractions as well as radial temperature profiles were published for the n-decane experiment. Figure 10a shows the soot volume fraction along the burner axis. It can be seen that the soot volume fraction maximum is underestimated by a factor of 1.6, about $6\,\hbox {mm}$ downstream of the experimental maximum. Deviations occur in the radial trends and values of the soot volume fraction for different heights above the burner as well, as can be seen in Fig. 10b. Profiles close to the burner entrance ($h=30\,\hbox {mm}$, $40\,\hbox {mm}$ and $50\,\hbox {mm}$) seem to reflect the radial trend very well, but noticeably overestimate the experimental measurements. The source terms for inception, growth and $\hbox {PAH}$ collision coincide at the respective maximum positions (compare Fig. 9c), leading to the excess of soot formation. The offset of the axial trend visible in Fig. 10a is also reflected in the radial trends for a height of $60\,\hbox {mm}$ and $70\,\hbox {mm}$. Although the overall soot yield is predicted considerably better, the trend along the radius diverges from the experiment. Due to the overpredicted soot volume fraction up to a height of $50\,\hbox {mm}$ associated radiative heat losses are also too high and the corresponding temperatures are slightly below the measurements, as can be seen in Fig. 10c. As soot volume fraction is increasingly underpredicted further downstream, the temperature is reproduced better at a height of $60\,\hbox {mm}$ and $70\,\hbox {mm}$. In the numerical study of Saffaripour et al. (2014), similar radial and axial trends of the soot volume fraction are predicted, even with a much more detailed gas phase and soot model. Measurement uncertainties in the experiment and not well predicted $\hbox {PAH}$ concentrations were identified as a cause for these differences. Another possible hypothesis for the observed deviations is given in the studies of Kholghy et al. (2013) and Jiang et al. (2019). Concurrently to the formation of $\hbox {PAHs}$ a formation of $\hbox {aromatic-aliphatic linked hydrocarbons (AALHs)}$ is found in the works. These liquid-like particles turn into soot by carbonization at temperatures of about $1500\,\hbox {K}$. This results in less $\hbox {PAHs}$ and soot close to the burner inlet and more soot near the center axis further downstream. However, it has to be mentioned, that this is only one explanation and up to now no model efforts on $\hbox {AALH}$ pathways are known to the authors.

Fig. 11

Measured and calculated profiles in an atmospheric laminar kerosene diffusion flame: a axial ${\mathrm{C}_2\mathrm{H}_x}$ mole fraction, b axial $\hbox {CH}_{4}$, $\hbox {CO}$ and $\hbox {CO}_2$ mole fraction, and c radial soot volume fraction. Open symbols: experimental data (Saffaripour et al. 2011); solid lines with filled symbols: calculated data; dotted lines with filled symbols: calculated data with a different kerosene surrogate (posf)

Temperature measurements are not available for the kerosene flame. Instead, axial species profiles of combustion intermediates ($\hbox {CH}_{4}$, $\hbox {C}_{2}\hbox {H}_{2}$, $\hbox {C}_{2}\hbox {H}_{4}$, $\hbox {C}_{2}\hbox {H}_{6}$) and products ($\hbox {CO}$, $\hbox {CO}_2$) are provided as shown in Fig. 11a, b. The calculated values are in an excellent agreement with the experiment. It demonstrates a good performance of the developed gas phase mechanism, which is essential for the quality of the $\hbox {PAH}$ and soot model, as these are strongly depending for example on $\hbox {C}_{2}\hbox {H}_{2}$ and $\hbox {OH}$. Figure 11c shows measured and simulated radial soot volume fraction profiles. At a height of $31\,\hbox {mm}$ an excellent agreement is obtained. In contrast to the n-decane flame the soot volume fractions above this height are significantly stronger underestimated. The chosen amount of aromatic components in the fuel could be critical for this. Radial trends agree well with increasing height. In addition to the deviations arising from the soot model, there also remains an influence from the investigated Jet A-1 mixture which is specified by the experiment with $30\%$ unidentifiable components. To investigate the impact of the chosen kerosene surrogate, the simulation is repeated with a different ÿÿcomposition based on POSF4658 (${39.7}\,\hbox {vol-}\%$ n-decane, ${30.7}\,\hbox {vol-}\%$ iso-octane, ${26.6}\,\hbox {vol-}\%$ toluene and ${3}\,\hbox {vol-}\%$ $\hbox {PAH}_{1}$) (Dooley et al. 2010). The higher content of aromatic components in this surrogate increases the soot volume fraction by a factor of 2 in the regions of the flame wings (see dotted lines in Fig. 11c) and could thereby explain the underestimation of the soot volume fraction peak at a height of $41\,\hbox {mm}$. However, regions close to the axis remain nearly unchanged and thus, the surrogate does not improve the results significantly.

As pointed out in the n-decane flame, the deviations between measured and predicted soot volume fractions on the burner axis may be related to uncertainties in boundary conditions and measurements, or result from missing $\hbox {AALH}$ pathways in the soot model. Compared to the numerically more expensive studies of Saffaripour et al. (2011) with 304 species and 2256 reactions, the presented model obtains slightly better results in soot volume fraction and comparable ones in the species profiles, for this test case. Even with the clearly more detailed model in their work, the soot volume fraction on the axis is not correctly reproduced. Taking accuracy and computational effort into account, the model performs extremely well over a wide range of fuels and operating conditions and thus, is well suited to be used in aero engine simulations.

4 Summary and Conclusions

A gas phase, $\hbox {PAH}$ and soot model for complex fuels has been presented for application in $\hbox {CFD}$ simulations. The developed model components, as well as significant results are:

1. The gas phase is described by an optimized kinetic scheme including 112 species and 846 reactions. It allows the simulation of short- and long-chain hydrocarbons and even fuel blends, e.g., kerosene. The mechanism achieves an excellent agreement with measurements of ignition delay times and species profiles. Even for fuel mixtures such as kerosene, the mechanism and the chosen three-component surrogate reflect important $\hbox {PAH}$ and soot model intermediates accurately. Compared to much more detailed kinetic schemes (several hundred species and many thousands of reactions) this gas phase offers a good tradeoff between the level of detail and computational effort.

2. An efficient sectional approach is used to model the $\hbox {PAHs}$. The reversibility of the $\hbox {PAH}$ gas phase interaction is taken into account by considering $\hbox {PAH}$ radicals. Soot is modeled with a sectional approach as well. The processes of inception, growth, $\hbox {PAH}$-soot collision, agglomeration and aggregation are covered. The good predicted sooting tendencies of various fuels in shock tube experiments demonstrate the flexible nature of the model. Especially trends at varying temperatures, pressures, and initial compositions are reproduced very well. The impact of polyynes on soot in the high temperature range has been investigated. It may provide a further path to soot formation at aero engine operating conditions, but the results obtained in this work do not allow a clear statement.

3. The smoking and non-smoking characteristics of an atmospheric laminar ethylene diffusion flame are well reflected. In addition, the model achieves good results for n-decane and kerosene flames. Soot intermediates, as well as oxidation species are accurately predicted. However, there is an underestimation of soot on the burner axis, which has to be addressed in the future. Additionally, the influence of premixing on soot formation is examined in Online Resource 1.