Phenalenyl growth reactions and implications for prenucleation chemistry of aromatics in flames

Abstract

The energetics and kinetics of phenalene and phenalenyl growth reactions were studied theoretically. Rate constants of phenalene and phenalenyl H-abstraction and C₂H₂ addition to the formed radicals were evaluated through quantum-chemical and rate-theory calculations. The obtained values, assigned to all π radicals, were tested in deterministic and kinetic Monte Carlo simulations of aromatics growth under conditions of laminar premixed flames. Kekulé and non-Kekulé structures of the polycyclic aromatic hydrocarbons (PAHs) evolving in the stochastic simulations were identified by on-the-fly constrained optimization. The numerical results demonstrated an increased PAH growth and qualitatively reproduced experimental observations of Homann and co-workers of non-decaying PAH concentrations with nearly equal abundances of even and odd carbon-atom PAHs. The analysis revealed that the PAH growth proceeds via alternating and sterically diverse acetylene and methyl HACA additions. The rapid and diverse spreading in the PAH population supports a nucleation model as PAH dimerization, assisted by the non-equilibrium phenomena, forming planar aromatics first and then transitioning to the PAH–PAH stacking with size.

---

Introduction

The formation of carbonaceous particulates in gaseous environments is of practical importance, a scientific curiosity, and a theoretical challenge. The gas-phase production of carbon black¹ and nanotubes² is currently at the industrial scale and manufacturing of graphene could become a reality.^3,4 Large-scale production of hydrogen coupled to carbon formation in plasma decomposition of methane is envisaged as a viable approach to a clean energy economy.⁵ The mechanism of particulate formation occurring in terrestrial gas-phase environments on a millisecond timescale has become a subject of research since Michael Faraday's brilliant demonstration of the phenomenon with a candle flame.⁶ Similar mechanistic characteristics were suggested to underlie interstellar dust formation.⁷ The quest for mechanistic details of the phenomena of carbon particle formation challenges the quantum-chemical and reaction-rate theories.^8,9 The present study is concerned with soot formation in hydrocarbon combustion.

Many aspects of soot formation are generally understood in mechanistic terms.¹⁰ One notable exception is the specifics of soot particle inception. The topic was recently reviewed,^11–13 and the present understanding can be summarized as follows. Polycyclic aromatic hydrocarbons (PAHs) are molecular precursors. They form from the initial fuel and their molecular fragments and grow afterwards through chemical reactions taking place at their edges. At some point, PAH moieties start forming clusters. The remaining principal questions aim at identification of the PAH (size, state, etc.) that initiates the clustering and its underlying mechanism. Thus, identification of PAHs and their molecular evolution preceding and surrounding the initial phase of clustering becomes imperative.

Our recent study examined the evolution of PAH molecular structures starting with naphthalene in detailed, sterically resolved kinetic Monte-Carlo (kMC) simulations.¹⁴ The numerical analysis was performed for the environments of two premixed flames of ethylene and the reaction set included results of recent PAH reaction studies. One of the major outcomes of the simulations was the prominence of phenalene (1H-phenalene) and its radical, phenalenyl. The present study expands on this observation and its implications for PAH evolution.

Phenalenyl is an odd-carbon-atom π-conjugated radical that has drawn scientific attention due to its high molecular symmetry and high thermodynamic stability.^15,16 The phenalenyl radical and its parent molecule, 1H-phenalene, have been observed in flames^17–19 and their formation reactions were explored theoretically.^20,21 Of interest to soot formation in flames is the phenalenyl propensity to dimerize,^15,22 as a possible initiating step of particle inception,^18,19 and its further growth into larger PAH structures enabling such inception.^14,23

Our prior kinetic study¹⁴ included growth of phenalenyl but in a rather generic form, treating reactions of phenalenyl-type structures as those of “regular” PAHs; in other words, all PAH-edge growth reactions were assigned to the same reaction-rate classes. Given the abundance of phenalenyl revealed in the prior numerical observations, here we examine these reactions with further scrutiny. First, we evaluated the rate constants for the following phenalene and phenalenyl reactions:

A2AH + H ⇆ πA2A + H₂ (1)

πA2A + C₂H₂ ⇆ A2AR5 + H (2)

πA2A + H ⇆ πA2A˙ + H₂ (3)

πA2A˙ + C₂H₂ ⇆ πA2AR5 + H (4)

where A2AH designates 1H-phenalene, πA2A is the phenalenyl radical, πA2A˙ is the π–σ phenalenyl biradical, and A2AR5 and πA2AR5 are respective cyclopenta adducts (see the identification of molecular structures in potential energy diagrams shown in Fig. 1 and 2). Reaction (1) was recently examined experimentally and theoretically at low temperatures;²⁴ here, we supplement these data for higher temperatures pertinent to combustion.

Fig. 1 Potential energy diagram for the C₁₃H₉ (phenalenyl, πA2A) + C₂H₂ reaction. All relative energies given in kcal mol⁻¹ with respect to the initial reactants are calculated at the G3(MP2,CC)//B3LYP/6-311G(d,p)+ZPE(B3LYP/6-311G(d,p)) level of theory.

Fig. 2 Potential energy diagram for the C₁₃H₈ (πA2A˙) + C₂H₂ reaction. All relative energies given in kcal mol⁻¹ with respect to the initial reactants are calculated at the G3(MP2,CC)//B3LYP/6-311G(d,p)+ZPE(B3LYP/6-311G(d,p)) level of theory. Black and blue curves correspond to the triplet and singlet electronic states, respectively.

The computed rate constants were then used in deterministic and stochastic simulations, performed for the same two flames examined earlier.¹⁴ The deterministic simulations were performed by solving differential equations for a naphthalene-to-pyrene reaction subsystem that allowed us to numerically compare growth via the phenalenyl and phenanthrenyl radicals. A broader analysis was performed by carrying out stochastic simulations that employed the kMC model of our prior study.¹⁴ The computed rate-constant values of reactions (1)–(4) were assigned, on the per-site basis, to the reaction classes of all phenalenyl-type sites appearing on PAH edges. The obtained numerical results confirm and expand on the conclusions of the prior study¹⁴ regarding the role of phenalenyl-type chemistry. Among other things, the kMC results reveal a pattern of PAH evolution consistent with experimental observations of Homann and co-workers,^25,26 yet offer a differing explanation for the underlying mechanism of the observed phenomena. A broader discussion on the implication of the present results to the subject of soot particle nucleation is given after the following description of the calculation methods and presentation of the numerical results.

Calculation methods

Quantum chemistry

To explore potential energy surfaces (PESs) of H abstractions from phenalene and phenalenyl radicals and C₂H₂ additions to πA2A and πA2A˙, we exploited our standard theoretical approach to electronic structure calculations used in numerous previous studies of various reactions involved in the PAH growth. In particular, density functional theory (DFT) B3LYP^27,28 calculations with the 6-311G(d,p) basis set were carried out for geometry optimization of all local minima and transition states and for computation of vibrational frequencies for the stationary structures. The B3LYP/6-311G(d,p) computed frequencies were utilized for the evaluation of zero-point vibrational energy (ZPE) corrections and in rate constant calculations. Single-point energies were refined using the composite model chemistry G3(MP2,CC) scheme^29,30 where the CCSD(T)/6-31G(d) energy is upgraded by adding a basis set correction obtained by MP2 calculations with the G3Large and 6-31G(d) basis sets and ZPE is also included. The B3LYP calculations were performed with the Gaussian 16 program package³¹ and the CCSD(T) and MP2 calculations were carried out with the MOLPRO 2021 code.³²

Minimal energy singlet–triplet crossing points (MSX) were initially optimized employing the MOLPRO package using the multireference complete active space CASSCF/6-311G** method^33,34 with the active space including ten electrons on 10 orbitals (10,10). Next, the MSX structures were re-optimized at the DFT ωB97XD/6-311G** level³⁵ using the NST code.³⁶ It should be noted that the optimized CASSCF and ωB97XD MSX geometries appeared to be very similar. NST, in addition to the MSX optimization, computes the Landau–Zener nonadiabatic transition state theory (TST) reactive flux for a transition between electronic states of different multiplicities via the MSX.^37–39 The spin–orbit coupling constants required for the NST calculations were computed at the same CASSCF(10,10)/6-311G** level as the geometry optimization of the MSX using MOLPRO. Single-point energies of the optimized MSX structures were refined by G3(MP2,CC) calculations of the triplet states.

Reaction rate constants

Temperature- and pressure-dependent rate constants in the temperature range of 500–2500 K and in the pressure range of 0.01–100 atm were computed employing the RRKM-ME approach as implemented in the MESS code,^40,41 generally within the rigid rotor-harmonic oscillator approximation (RRHO). Internal rotations corresponding to low-frequency vibrational modes were treated as hindered rotors and their internal rotational potentials were evaluated using B3LYP/6-311G(d,p) PES scans along the corresponding torsional angles. The “exponential down” model⁴² was used to treat the collisional energy transfer in ME, with the temperature dependence of the range parameter α for the deactivating wing of the energy transfer function expressed as α(T) = α₃₀₀(T/300 K)ⁿ. The values of n = 0.62 and α₃₀₀ = 424 cm⁻¹ in this expression were earlier derived from classical trajectory calculations^43,44 and, along with the Lennard-Jones parameters ε = 390 cm⁻¹ and σ = 4.46 Å, were adopted from our earlier studies of the growth kinetics of two- and three-ring PAHs, including the studies of C₂H₂ addition reactions to various PAH radicals.^21,45–50

Reaction kinetics

The computed rate constants were tested in two sets of calculations, deterministic (ODE) and stochastic (kMC). Both sets of simulations were performed in the environments of the two atmospheric burner-stabilized flames of ethylene considered in our prior study,¹⁴ namely a stagnation 16.3% C₂H₄-23.7% O₂–Ar flame of Wang and co-workers⁵¹ (cold gas velocity 8.0 cm s⁻¹ and burner-to-stagnation surface separation 0.8 cm; designated as C3H08) and a 15.6% C₂H₄-17.7% O₂–N₂ flame of Faeth and co-workers⁵² (designated as XSF1.88). These are well-studied and well-understood flames, ranging from low (C3H08) to heavier (XSF1.88) soot loading. The two flames were simulated with the FFCM1 reaction model⁵³ using the Cantera software.⁵⁴ The computed flame temperature and species concentration profiles of H, H₂, CH₃, C₂H₂, O, OH, and O₂ were supplied to the ODE and kMC solvers.

The deterministic simulations were performed by solving a set of ordinary differential equations (ODE) describing a reaction submodel of naphthalene growth into pyrene via two primary competing reaction pathways, those via phenalenyl and phananthrenyl radicals. The ODE system was solved using the Matlab ode15s stiff integrator.⁵⁵ The initial concentration of naphthalene in these calculations was assigned a value computed in the flame simulations described above using Appel et al.'s (ABF) reaction model.⁵⁶ Qualitatively similar results, for the purpose of the present study, were obtained using other reaction models. Unfortunately, the prediction of naphthalene varies substantially among these models, reflecting the present uncertainty in PAH chemistry in general and in naphthalene formation reactions specifically (see, e.g., Valencia-López et al.⁵⁷).

The kMC simulations tracked a single PAH molecular structure evolving in the flame environment. Each kMC run started with a naphthalene (A2) molecule “placed” in a flame location at 1400 K and carried out for a duration of 1 ms. The stochastic evolution of the PAH structure was simulated with the previously developed Matlab code¹⁴ using the Gillespie algorithm.^58,59 The first-order per-site reaction rate constants were calculated using the time-dependent temperature and gaseous species profiles obtained in the flame simulations. The rate-constant values were updated every 10 μs. Two sets of simulations were performed, 100 000 runs for the C3H08 flame and 20 000 runs for the XSF1.88 flame, thus obtaining for each case about 1000 runs surviving early destruction of the aromatic substrate through oxidation and thermal decomposition.

Every PAH structure emerging in the present kMC simulations was identified as being Kekulé or non-Kekulé and the assignment of the pertinent reaction event and its rate was based on this instantaneous determination. The latter was accomplished through constrained optimization, as follows. The evolving aromatic structure in the kMC code¹⁴ is represented by undirected graph G(V,E) with n vertices V representing the carbon atoms and m edges E representing the C–C bonds.⁶⁰ We introduce vector whose elements are edge multiplicities, equal to 1 for a single C–C bond and 2 for a double C=C bond. The vertex valency is the sum of its incident-edge multiplicities. The kekulization problem can then be stated as minimization of the number of vertices whose valency deviates from the “Kekulé” valency, i.e., 4 for the inner vertices (with vertex valency 4) and 3 for the outer vertices (with vertex valency 3). Mathematically, this takes the form

where A is an n × m matrix with elements A_ij = 1 for edge j being an incident edge of node i and A_ij = 0 otherwise,and b is an m × 1 column vector whose elements are “Kekulé” valences of the corresponding vertices. The optimum x gives the n × 1 vector, b − Ax_opt, whose non-zero elements (equal 1) indicate vertices (C atoms) with the “non-Kekulé” vacancy. Thus, all elements being zero indicates a Kekulé structure and non-Kekulé otherwise. The optimization was solved by Matlab mixed-integer linear programming function intlinprog. It was found to execute sufficiently fast for PAH structures up to about 100 carbon atoms. Above this size, the criterion for the Kekulé/non-Kekulé structure was the even/odd number of carbon atoms, possibly introducing a very small error to simulations of the larger structures, based on the results discussed later in the text.

Results and analysis

Potential energy surfaces, rate constants and equilibrium constants

Let us first consider reactions (1) and (2), the hydrogen-abstraction-C₂H₂-addition (HACA) growth¹⁰ with the phenalene molecule (A2AH) was tested. In reaction (1), H abstraction from the CH₂ group in phenalene producing the phenalenyl radical πA2A is computed to be highly exothermic, by 41.4 kcal mol⁻¹, and to feature a low barrier of 2.7 kcal mol⁻¹. These values are close to the earlier predictions from CCSD(T)/cc-pVTZ//B3LYP/cc-pVTZ calculations by Schneiker et al.,²⁴ 38.7 and 2.5 kcal mol⁻¹, respectively. The high exothermicity and low barrier of the reaction originate from a very weak C–H bond in A2AH. Fig. 3(a) illustrates the computed rate constants of reaction (1) in the forward and reverse directions. The forward reaction appears to be orders of magnitude faster than the reverse one reaching ∼1 × 10⁻¹⁰ cm³ molecule⁻¹ s⁻¹ at 1500 K. The equilibrium constant is 7.8 × 10⁵ at 1500 K, reducing to 2.5 × 10³ at 2500 K. Under the flame conditions simulated in the present study, the reduced equilibrium constant of this reaction, which is the product of the equilibrium constant and the H-to-H₂ concentration ratio,¹¹ is greater than 70 and in the soot nucleation part of the flame it is greater than 100. These values indicate that reaction (1) is practically irreversible under the conditions of interest and A2AH is replaced by πA2A.

Fig. 3 Calculated TST and RRKM-ME rate constants for: (a) the direct H abstraction reactions A2AH + H ⇆ πA2A + H₂ (1) and πA2A + H ⇆ πA2A˙ + H₂ (3) in the forward and reverse directions; (b) the πA2A + C₂H₂ ⇆ A2AR5 + H reaction (2) including the total and individual channels at 1 atm; (c) the πA2A˙ + C₂H₂ ⇆ πA2AR5 + H reaction (4) including the total and individual channels at 1 atm; (d) the reverse p1′′ + H reaction including the total and individual channels at 1 atm.

The second step of the HACA growth sequence is the addition of C₂H₂ to the formed radical site followed by an extra ring closure and possibly an H atom loss restoring aromatization. However, such reaction (2) appears to be unfavourable for πA2A (Fig. 1). Phenalenyl, as any other PAH π-radical is unreactive and the acetylene addition step to form a complex i1′ features a high barrier of 27.9 kcal mol⁻¹. Next, i1′ can either lose a hydrogen atom to form ethynyl-substituted phenalene p1′ or undergo cyclization. The H loss channel to p1′ is overall endothermic by 45.0 kcal mol⁻¹, with the critical barrier at a H elimination step of 50.4 kcal mol⁻¹ relative to the πA2A + C₂H₂ reactants. The cyclization may occur immediately in i1′ forming i2′ or be preceded by H migration i1′ → i4′ eventually leading to i3′. Both i2′ and i3′ include an additional five-membered ring added to the phenalene core (A2AR5H), are connected through a 1,2-H shift within this five-membered ring, and can split a hydrogen atom to produce two isomers of A2AR5, p2′ from i2′ and p3′ from i3′. Here, the i1′ → ts3′ → i2′ → ts4′ → p2′ + H pathway with the highest in energy transition state at 24.6 kcal mol⁻¹ above πA2A + C₂H₂ (ts4′) is more advantageous than either i1′ → ts3′ → i2′ → ts5′ → i3′ → ts6′ → p3′ + H or i1′ → ts7′ → i4′ → ts8′ → i3′ → ts6′ → p3′ + H. Thus, although the p3′ + H products are thermodynamically favourable by 14.1 kcal mol⁻¹, the p2′ + H products are expected to be preferable kinetically. This is indeed confirmed by the rate constant calculations, as can be seen in Fig. 3(b). In the 1000–2250 K range, reaction (2) is predicted to form p2′ + H, whereas at higher temperatures the production of ethynyl-substituted phenalene p1′ + H becomes slightly more favourable. However, due to the high barriers involved, reaction (2) is very slow, with a calculated rate constant of only 2 × 10⁻¹⁷ cm³ molecule⁻¹ s⁻¹ at 1500 K and 1 atm. This value is nearly five orders of magnitude lower than, for example, the rate constant for C₂H₂ addition to a σ-PAH radical phenanthrenyl leading to the formation of pyrene + H. We discussed a slow reactivity of resonantly stabilized π-PAH radicals toward acetylene in a previous work⁶¹ and here phenalenyl πA2A follows the predicted behaviour. Thus, the direct PAH enlargement involving πA2A via C₂H₂ addition is unlikely to be competitive.

However, can πA2A undergo the HACA growth through a different path? Since, according to the calculated equilibrium constant for reaction (1), phenalene mostly exists in sooting flames in the form of πA2A, it can partake in other bimolecular reactions and in particular be subjected to a direct hydrogen abstraction reaction with an H atom, reaction (3). In this case, H abstraction is preferable from one of six equivalent positions on the edges of πA2A and produces a π–σ biradical πA2A˙ with a triplet ground electronic state and one unpaired electron located on a π-orbital and the other on a σ-orbital. Reaction (3) is predicted to be 7.8 kcal mol⁻¹ endothermic with a barrier of 15.5 kcal mol⁻¹ in the forward direction. This energetics is comparable with endothermicity and the barrier height for the prototype C₆H₆ + H → C₆H₅ + H₂ reaction activating benzene to the phenyl radical, 8.8 and 17.0 kcal mol⁻¹, respectively, at the same G3(MP2,CC)//B3LYP level of theory.⁶² The forward and reverse rate constants for reaction (3) are depicted in Fig. 3(a). The forward direction is preferable above 1000 K, with the equilibrium constant reaching 3.4 at 1500 K and increasing to 8.1 at 2500 K. The calculated rate constant is comparable with that for C₆H₆ + H, which has the same reaction path degeneracy as reaction (3). In particular, at 1500 K, the calculated rate constants are 3.3 × 10⁻¹² cm³ molecule⁻¹ s⁻¹ for reaction (3) and 4.9 × 10⁻¹² cm³ molecule⁻¹ s⁻¹ for C₆H₆ + H.⁶³ A conclusion thus can be made that PAH π-radicals are expected to behave similarly to closed-shell PAH molecules in a sense that they are generally unreactive, have long lifetimes, and can be activated through direct H abstraction reactions with rate constants that are alike. The difference is that H abstraction from closed-shell PAH molecules produces σ-PAH radicals, whereas that from the π-radicals generates π–σ biradicals like πA2A˙.

The next step in the HACA mechanism is the acetylene addition to πA2A˙ initiated on a triplet potential energy surface (Fig. 2). Here, the barrier for the C₂H₂ addition to the σ-radical site in the ring is low, only 2.5 kcal mol⁻¹ and the initial complex produced, i1′′t, is stabilized by 40.9 kcal mol⁻¹ relative to the πA2A˙ + C₂H₂ reactants. Several reaction scenarios are possible commencing from the triplet isomer i1′′t. First, it can lose the hydrogen atom from the α-carbon in the side chain leading to the ethynyl-substituted phenalenyl radical p1′′. The H loss proceeds via ts9′′ located 4.0 kcal mol⁻¹ below the reactants and the p1′′ + H products are 19.7 kcal mol⁻¹ exothermic. Second, the ring closure pathway on the triplet PES leads to i2′′tvia ts3′′t (−15.4 kcal mol⁻¹) and i2′′t can either directly eliminate a hydrogen atom forming a resonantly stabilized π-radical πA2AR5 (p2′′) or first undergo a 1,2-H migration to i4′′t followed by an H loss to p2′′. Also, i1′′t can feature a 1,5-H shift to i3′′t prior to the ring closure to i4′′t. Among the two pathways leading to πA2AR5 on the triplet surface, i1′′t → ts5′′t → i3′′t → ts6′′t → i4′′t → ts7′′t → p2′′ + H is more energetically preferable with the highest barrier at ts6′′t (−25.6 kcal mol⁻¹ relative to the reactants) than either i1′′t → ts3′′t → i2′′t → ts4′′t → p2′′ + H or i1′′t → ts3′′t → i2′′t → ts8′′t → i4′′t → ts7′′t → p2′′ + H both proceeding via the critical ts3′′t. In addition to the triplet PES, there are also reaction routes involving a singlet surface after triplet–singlet intersystem crossing (ISC) takes place. We have located two minima of the seam of crossing, one in the vicinity of i1′′t and the other near i3′′t. Isomer i1′′ includes a carbene CH group and can exist both in the triplet (−40.9 kcal mol⁻¹) and in the singlet (−38.9 kcal mol⁻¹) electronic states. The msx1 structure features a similar geometry and its energies are 3.8 and 1.8 kcal mol⁻¹ higher than those of i1′′t and i1′′s, respectively. The singlet carbene structure formed after ISC undergoes facile ring closure to i2′′svia a low barrier of 4.6 kcal mol⁻¹ (ts3′′s, −34.6 relative to the reactants). i2′′s represents one of the A2AR5 isomers and a potential product (p2′) of reaction (2). Isomer i3′′ has a lone pair on one of the ring carbon atoms and also has close in the energy triplet (−40.2 kcal mol⁻¹) and singlet (−37.7 kcal mol⁻¹) states. The ISC process takes place at msx2 which lies at −35.0 kcal mol⁻¹ relative to the reactants. Once i3′′s is formed after ISC, it should easily ring-close to i4′′s, the other A2AR5 isomer p1′. The singlet isomers i1′′s and i4′′s lie much lower in energy than any of the structures on the triplet surface and thus the ISC process is overall highly exothermic. However, as will be discussed below, the crossing to the singlet surface does not make a significant contribution to the reaction kinetics due to low spin–orbit coupling constants at msx1 and msx2. Both A2AR5 structures i1′′s and i4′′s can lose a hydrogen atom producing πA2AR5 p2′′ without exiting barriers. The energies of the weakest C–H bonds in i1′′s and i4′′s are rather low, 52.0 and 66.1 kcal mol⁻¹, respectively, which is comparable with the C–H bond strength in the CH₂ of phenalene, 62.2 kcal mol⁻¹. Therefore, similarly to phenalene, both A2AR5 isomers under flame conditions, i.e., at high temperatures and with significant concentrations of H atoms and other H-abstracting radicals, are likely to irreversibly convert to the π-radical πA2AR5. Summarizing, the πA2A˙ biradical can grow an extra five-membered ring on its zig-zag edge by adding acetylene thus propagating the chain of π-radicals from πA2A to πA2AR5 via the HACA mechanism. It is pertinent to note that the 165 and 189 m/z signals have been commonly seen in mass spectra in numerous flame studies^64–66 and, in particular, they were hypothesized to correspond to πA2A to πA2AR5.⁶⁷

The rate constants for reaction (4) calculated at 1 atm are depicted in Fig. 3(c) (a complete set of the computed rate coefficients is given in the ESI†). The reaction is relatively fast, with the rate constant ranging between ∼10⁻¹² and 10⁻¹¹ cm³ molecule⁻¹ s⁻¹ in the 750–2500 K temperature interval. These values are comparable with the rate constants for acetylene addition to PAH σ-radicals.⁶¹ At low temperatures, the preferable reaction product is collisionally stabilized i4′′t, but at the temperatures of interest to combustion, the reaction is dominated by the formation of πA2AR5 p2′′. This product is favourable above 800 K and up to ∼2000 K, whereas at the highest temperatures, the yield of the ethynyl-substituted phenalenyl radical p1′′ slightly exceeds that of p2′′. Nevertheless, in the reverse p1′′ + H reaction, p1′′ mostly undergoes H-assisted isomerization to p2′′ (Fig. 3(d)). In this reaction, collisional stabilization of i4′′t again prevails at low temperatures, p2′′ + H exhibits the highest branching ratio at 900–2000 K, whereas the return to the πA2A˙ + C₂H₂ reactants takes over as the preferred channel above 2000 K. These results corroborate the conclusion that πA2A can efficiently undergo a HACA growth by adding an extra five-membered ring on a zig-zag edge.

It is interesting to consider the role of ISC and the singlet PES in reaction (4). This role is minor as indicated by the fact that the rate constant for collisional stabilization of the much more favourable singlet isomer i4′′s is up to an order of magnitude lower than that for its triplet counterpart i4′′t. More particularly, Fig. 4 compares high-pressure limit rate constants for ISC of i1′′t and i3′′t with those for their isomerization reactions occurring on the triplet surface without changing the spin multiplicity. Despite the fact that “activation energies” for ISC, 3.8 and 5.2 kcal mol⁻¹ for msx1 and msx2, respectively, are significantly lower than the isomerization barrier heights, 9.6/8.9 kcal mol⁻¹ for i1′′t ⇆ i3′′t and 14.6 kcal mol⁻¹ for i3′′t → i4′′t, the ISC rate constants are normally few orders of magnitude lower than their isomerization counterparts. This is caused by the fact that the computed spin–orbit coupling constants are as low as 1.5 and 3.0 cm⁻¹ for msx1 and msx2, respectively. The only range of temperatures at which the i3′′t → i3′′s step is competitive with i3′′t → i4′′t is below 600 K. While the isomerization rate constants show a typical Arrhenius behaviour steadily growing with temperature, the ISC rate constants slightly increase up to 1250–1400 K and higher temperature calculations of the ISC rates require a refinement of the energy grid and were not pursued because of the insignificance of these channels. It is clear that in the relevant temperature interval, the reaction predominantly proceeds on the triplet PES without crossing onto the singlet surface.

Fig. 4 RRKM-calculated rate constants for isomerization and intersystem crossing of the i1′′t and i3′′t intermediates in the high-pressure limit.

Reaction kinetics

The deterministic (ODE) calculations were performed with the reaction model composed of 41 reactions, describing the transformation of naphthalene into pyrene. Only reaction paths (3) and (4) were included; reaction (1) is essentially shifted to the right, as discussed above, and numerical tests confirmed no meaningful flux was observed by including reaction paths (1) and (2). The principle reaction pathways are illustrated in Fig. 5 and the complete set of reactions and their rate constants are presented in Table S2 of the ESI.†

Fig. 5 Principal reaction pathways tested in the ODE simulation: phenalenyl is marked with red arrows and phenanthrene is marked with blue arrows.

We refer to the reaction pathways shown in Fig. 5 as principle because one of them, marked with blue arrows, is the HACA growth with acetylene being the growth species and includes only even carbon-atom intermediates, whereas the second pathway, marked with red arrows, is the HACA growth with acetylene and methyl being the growth species and includes even and odd carbon-atom intermediates. The numerical results presented in the middle panels of Fig. 6 show that the additional growth pathways via π-radicals increases the production of pyrene, as compared to the previous model of only even carbon-atom growth.⁴⁹ Also, the concentrations of the pathway intermediates, especially of pyracyclene and cyclopenta-phenalenyl, are significantly higher than that of pyrene, as demonstrated in the bottom panels of Fig. 6. These results, consistent with our prior numerical observations,¹⁴ suggest growth proceeding via multitude of sterically enabled pathways. The latter was explored in the kMC simulations described next.

Fig. 6 Left panels: C3H08 flame; right panels: XSF1.88 flame; top-row panels: flame temperature; middle-row panels: concentration of pyrene (A4) obtained with the prior mechanism⁴⁹ (designated as “base”, blue) and with the new reactions added (red); bottom-row panels: computed concentrations of pyracyclene (R5A2R5, black), cyclopenta-phenalenyl (πA2AR5, red), phenalenyl (πA2A, green), and pyrene (A4, blue).

The kMC simulations were carried out using the reaction model of the previous study¹⁴ with two additional reaction classes, chemically similar to reactions (3) and (4) and with their computed rate constant values assigned, on the per-site basis, to all reaction events of these classes. Identification of Kekulé structures, implemented in the present kMC code, allowed making appropriate rate constant assignments.

Fig. 7 displays frequencies of PAH structures observed in the kMC simulations. The maxima of the PAH-structure counts for every number of carbon atoms are shown in the top panels, obtained by taking the maxima of PAH appearances over the entire flame simulation. This form of display was chosen to make a comparison to remarkable experimental observations of Homann and co-workers (displayed in Fig. 5 of Keller et al.²⁶) for the maximum concentrations of PAH collected over a laminar flame: (a) the odd and even carbon-atom concentrations are close to each other, (b) both exhibit a non-decaying pattern for sizes over 20–25 carbon atoms, and (c) both decay sharply for the smaller sizes. The results of the present kMC simulations, displayed in the top panels of Fig. 7, reproduce (qualitatively) all these experimental features.

Fig. 7 Computed frequencies of PAH structures displayed by the number of C atoms and reaction times in the two flames, C3H08 and XSF1.88.

The general agreement of the kMC results produced for the two substantively different ethylene–air flames studied in the present study¹⁴ and the concurrence of the computed patterns with the experimental observations in a very different flame (low-pressure benzene–oxygen) suggest the universality of these features.

Further confirmation of the phenomena inferred from the “summary” display, utilized in the top panels of Fig. 7 and in Fig. 5 of Keller et al.,²⁶ is presented in the middle and bottom panels of Fig. 7 by displaying, respectively, distribution of PAH sizes at several reaction times and abundances of several PAH sizes versus time. These results illustrate the development of the distributions of individual-mass PAH and demonstrate the consecutive appearance of their concentration peaks.

Discussion

Mechanism of odd and even carbon-atom aromatics growth

Homann and co-workers²⁶ suggested that their experimental observations provide evidence against the HACA mechanism promoting the PAH growth; our present results demonstrate the opposite.

The growth of PAH in high-temperature environments, starting with a single aromatic ring, faces thermodynamic resistance, i.e., the reversibility of carbon-adding reactions. The underlying feature of the HACA mechanism, discovered⁶⁸ through numerical simulations of PAH growth and explained^10,14 afterwards, is alleviation of the thermodynamic resistance through repetitive reactivation of aromatic-edge sites by hydrogen atoms preceding the addition of carbon carriers to these sites. Among possible carbon-atom carriers, those present in sufficient abundance and whose addition lowers the potential energy or increases the reaction entropy or both should dominate the second, carbon-addition step of the HACA sequence. In high-temperature hydrocarbon pyrolysis and combustion systems, such a species is primarily acetylene, especially when its addition releases a hydrogen atom. Obviously, other species, under particular conditions, could become competitive with acetylene. One such species is the methyl radical, CH₃.

The methyl radical was suggested as a possible PAH growth species as early as 1984 by Weissman and Benson.⁶⁹ However, its repetitive addition was shown to be not sufficiently rapid.⁷⁰ Yet, the addition of CH₃ expanding a five-member ring into a six-member ring has emerged as a viable process. The possibility of such a ring expansion was suggested independently from experimental observations,^71,72 computational tests,⁵⁹ and theoretical examination of the energetics and kinetics of elementary reactions.^21,45,73 Our preceding kMC study,¹⁴ employing the latest rate constants for the methyl-induced ring expansion, showed that the latter reactions are competitive with the acetylene-addition PAH growth. The results of the present study, employing a detailed pathway of odd carbon-atom growth, based on reactions (3) and (4), reaffirm the prior conclusion.

Both the methyl and acetylene growth reaction pathways are HACA, H-Abstraction-Carbon-Addition,^10,11 and hence will be referred to hereafter as HAC_H3A and HAC_2H2A, respectively. The kMC results indicate that the PAH growth proceeds as an alternation of HAC_2H2A and HAC_H3A steps, as illustrated in Fig. 8. Both even and odd carbon-atom aromatics retain their “parity” (even remain even or odd remain odd) undergoing the respective HAC_2H2A reaction steps, and change their parity (even turns into odd and odd into even) with the HAC_H3A steps. While the parity switching can possibly be accomplished with other reactions, the present results demonstrate that it is the parity switching that explains the Keller et al.²⁶ observation of nearly equal concentrations of even and odd carbon-atom PAHs. Also, the alternating HAC_2H2A and HAC_H3A sequence includes elements of prior suggestions.^45,59,72

Fig. 8 PAH growth through the alternation of HAC_2H2A and HAC_H3A reaction steps.

Homann and co-workers²⁶ interpreted their experimental observations as “throw[ing] doubt on the general validity of the so-called HACA mechanism,” yet the detailed numerical analysis of the present study accentuates the underlying reaction mechanism as HACA. Indeed, the present study is in accord with Keller et al.'s conclusion that “the delocalized unpaired electron does not increase the reactivity of a PAH π-radical over that of a closed-shell PAH with a similar structure.” In addition to the thermodynamic stability of the PAH π-radicals, like phenalenyl (πA2A), the theoretical results showed that the addition of H atoms, reaction (1), is shifted to the π-radicals and the additions of C₂H₂ molecules, reaction (2), is relatively slow. The growth of these π-radicals is attained via the HACA mechanism, HAC_2H2A or HAC_H3A: formation of a π–σ biradical through H abstraction, a lareaction (3), followed by carbon addition, a lareaction (4), both reactions having rates comparable to those of the PAH molecules.^46,49,63

Distribution of PAH concentrations

The present kMC results reproduce (qualitatively) the distribution of concentrations of aromatics observed by Kellers et al.²⁶ As can be seen in the top panels of Fig. 7, the initial PAH molecules and π-radicals, up to about the size of 25 C atoms, fall rapidly and then remain at about the same size up to about 120 C atoms. There is also a slight increase in size, similar to the one noted in the experiment.²⁶

The initial drop is a consequence (and thus the confirmation) of the HACA mechanism. Growth of the initial PAH is a “struggle” against reversibility of the reaction sequence. The first “break” in this struggle is the formation of acenaphthylene, the first island of stability.⁶⁸ The term “island of stability” signifies a chemical species that is not only thermodynamically stable but also kinetically stable, showing resistance to the destruction through thermal decomposition and oxidation of its radicals formed by H abstraction.

The next and a stronger island of stability was identified⁶⁸ to be pyrene; its formation is practically irreversible and thermal decomposition⁷⁴ and oxidation⁷⁵ of pyrenyl are relatively slow. Thus, the kinetic sequence of aromatic growth up to pyrene faces the largest struggle against the reversibility and, as a result, the sequential aromatic intermediates drop in the concentration, especially prior to acenaphthylene.

The present results indicate that phenalenyl is also an island of stability, as well as other phenalenyl-like π-radicals, assuming that π–σ biradical's decomposition and oxidation reactions are similar to those of pyrenyl. This implies that there is a path for growth past pyrene that moves continuously through islands of stability, as exemplified in Fig. 8. Furthermore, larger PAH structures offer a greater number of kinetically equivalent edge sites, thereby increasing the multiplicity of growth pathways. The increase in the multiplicity and diversity of closely spaced islands of stability evidently results in the non-decreasing concentrations of PAH observed by Kellers et al.²⁶ and reproduced by the present kMC simulations.

Implications to soot particle inception

The phenomenon of the non-decreasing concentrations of aromatic structures, starting at around 20–25 carbon-atom sizes, has significant implications to the currently debated subject of soot particle inception.

In the initial detailed modeling of soot formation,^56,76–78 particle nucleation was represented by PAH clustering: binary collisions of all PAHs past a prescribed size. The PAH growth from the prescribed size up to an infinite size was modeled by a replicating HAC_2H2A sequence^68,79 and its kinetics solved with a mathematically rigorous algorithm of linear lumping.^79,80 The PAH initiating the clustering was presumed to be (ace)pyrene^56,78,81 or coronene^76,77,82 and the collisions of all PAHs were assigned the efficiency of unity. It has been established since then that pyrene dimerization is thermodynamically unstable^83–86 yet starting the clustering with a larger PAH, whose dimer is attaining the stability, does not reproduce the timing of soot appearance in flames.¹¹ The model and alternative proposals were recently scrutinized,^11,12 bringing about two basic propositions.

One of these propositions is the formation of a doubly bonded PAH dimer via a reaction of a cyclopenta group of a PAH (like in acepyrene) with a radical edge site of another PAH, either at a zig-zag site¹¹ (like in pyrene) or at a bay site⁸⁷ (like in phenanthrene). The reaction product, termed as the E-bridge, is formed as an angled PAH dimer but rapidly flattens out through H abstraction.⁸⁸ The other proposition pursues the premise of PAH π-dimerization.^{12,18,23,89–93} It has been established,^{15,16,22,90,94} though, that π-monoradicals, such as phenalenyl, preferentially form edge σ-bonded dimers rather than π–π stacked ones, both bonds being relatively weak. As a resolution to this, Martin and co-workers^12,90,92 suggested dimerization of π-biradicals. Indeed, π-radicals are known to exhibit a polyradical character,^15,90,95 as illustrated by three molecular examples shown on the left of Fig. 9. However, only one instance of such a species, shown on the right of Fig. 9, was detected in the kMC simulations of the present study. Such a low frequency is unlikely to explain PAH dimerization in flames.

Fig. 9 PAH π-polyradicals. Left: Molecular structures of phenalenyl, C₁₃H₉, and its triangular extensions, C₂₂H₁₂ and C₃₃H₁₅, computed by the constrained optimization of the present study, eqn (5); right: the only π-biradical, C₄₈H₁₈, detected in the present kMC simulations.

These considerations suggest that dimerization of smaller-sized PAHs may contribute to the “planar” growth (cross-linking) of PAH and not to their immediate stacking. The extent of such contributions is yet to be determined, but even without them, the PAH growth in size is already rapid. Also, the chemical opportunities provided by the planar dimerization may end up having little consequences, as the formed planar dimers can undergo rapid rearrangements and overgrowth due to the ring migration and HAC_2H2A/HAC_H3A reactions.⁹⁶ In fact, the latter processes, on their own, can create molecular structures, as illustrated in Fig. 10, having five-membered-ring features attributed²³ to PAH–PAH cross-linking whereas no such PAH reactions were included in the model.

Fig. 10 Examples of PAH molecular structures formed in the kMC simulations of the present study.

Turning to the particle inception, we now represent the available knowledge by the following model. We assume that PAH–PAH dimerization starts with the first major island of stability (like pyrene or acepyrene) forming planar aromatics first and then transitioning to the PAH–PAH stacking with size. The PAH dimerization for all sizes is assisted by the non-equilibrium phenomenon⁹ and is treated as irreversible¹¹ with an assigned collision efficiency, γ. We represent the latter by a switch function

that transitions from zero to unity with the reduced mass of the colliding PAH pair, μ, expressed by their number of C atoms. Parameters a and b in eqn (7) are adjusted by fitting γ at two points. The first fitting point is γ(16 × 18/34) for a pyrene–acepyrene pair at 1500 K, computed as a collision efficiency of 1 × 10⁻³, determined for the pyrenyl–acepyrene dimerization,⁹ multiplied by the [pyrenyl]/[pyrene] ratio computed from the steady-state expression for the pyrenyl formation by the reversible H abstraction and H addition reaction systems;¹¹ for flame C3H08, γ_16,18 = 5.92 × 10⁻⁶ with a = 0.325, and for flame XSF1.88, γ_16,18 = 1.10 × 10⁻⁶ with a = 0.370. The second fitting point, γ(27) = 0.5 for both flames, is for a binary collision of the 54 carbon-atom circumcoronenes, presumed to be approaching the van der Waals stacking.^97–99

To numerically assess such a model, we calculated a ratio of the total PAH–PAH dimerization rate, summing up the binary collisions among all the produced PAHs with the corresponding γ values, to that of the pyrene–acepyrene collision rate,

where N_i is the kMC count of PAHs with i C atoms at a given time, μ_ij is the reduced mass of the i–j collision pair, and σ_ij is the collision cross-section, expressed by the number of C atoms, n,with the PAH diameter assumed proportional to the cubic root of the PAH mass.¹⁰⁰ Note that all constant multipliers in eqn (8) and (9) are cancelled out by taking the ratio, eqn (8). The computed values of are displayed in Fig. 11.

Fig. 11 PAH–PAH dimerization ratio, computed by eqn (8) for the two flames, C3H08 and XSF1.88. The black solid line is the flame temperature. The symbols are values computed with different starting PAHs in the summation of the numerator of eqn (8): blue circles—pyrene–acepyrene, green squares—coronene–acecoronene, and red pluses—pyrene–pyracyclene. In all the cases, the summation was carried out up to PAH with 300 carbon atoms, although the computed results essentially stopped changing after the summation was limited to 80 carbon atoms.

There are several features in these results to note. First and foremost is that the total dimerization rate rapidly increases above that of pyrene/acepyrene alone, indicating the importance of accounting for a larger pool of PAH than just pyrene/acepyrene. The assumption of pyrene being a sole nucleating species, unfortunately propagating though much of the literature, evidently has its origin in misapprehension of the source sited.^56,78,81 As mentioned above (and in the previous accounts^10,11), the above-cited initial model included dimerization of all PAHs after pyrene. The growth of PAHs in this initial model was described by the only known at the time HAC_2H2A reaction sequence, whereas the present study shows a substantively faster PAH growth via alternating and sterically diverse HAC_2H2A and HAC_H3A reactions.

A consequence of the rapid and diverse spreading in the PAH population is that nucleation does not occur with a particular species (like pyrene or circumcoronene), but it is a process diffused over a sequence of them. In fact, this is a beneficial feature of treating soot nucleation kinetically¹⁰¹ as compared to rather fruitless applications of the classical nucleation theory that focusses on the identification of the critical-size nucleus.¹⁰²

Three sets of computations that differ in the PAH size that initiates the dimerization count are displayed in Fig. 11. The total dimerization rates for all the three sets overlap with each other after some initial time, as the increasing-in-size PAH population begins to dominate the sum. During the initial period, a larger starting PAH, coronene, slows the dimerization, as was seen earlier,⁷⁸ and a smaller starting PAH, pyracyclene, enhances the dimerization, consistent with the previous study.¹⁴

While parts of the model need to be refined and new reactions are likely to be discovered, the presented results clearly show that the alternating and sterically diverse HAC_2H2A and HAC_H3A reaction mechanism is an integral and important part of soot particle nucleation.

Summary and conclusions

Rate constants of phenalene and phenylenyl reactions, including both H-abstraction and C₂H₂ addition to the radicals formed, were theoretically evaluated. Their inclusion in naphthalene-to-pyrene deterministic simulations resulted in a larger production of pyrene under the conditions of the two laminar premixed flames.

Assigning these rate constants to reaction classes of all π-radicals, stochastic simulations reproduced experimental observations of Homann and co-workers²⁶ of non-decaying PAH concentrations, for sizes over 20–25 C atoms, with nearly equal abundances of even and odd carbon-atom PAHs. However, in contrast to Homann et al.'s interpretation of their observations as providing evidence against HACA, our present study accentuates the underlying reaction mechanism as HACA.

The analysis of the kMC results revealed that the PAH growth proceeds via alternating and sterically diverse HAC_2H2A and HAC_H3A reactions. The rapid and diverse spreading in the PAH population provides the support to a nucleation model as PAH dimerization, assisted by the non-equilibrium phenomena,⁹ forming planar aromatics first and then transitioning to the PAH–PAH stacking with size.

Author contributions

AMM: quantum-chemical and reaction-rate calculations, conceptualization, analysis, and writing; AWJ: NST calculations; MF: ODE and kMC calculations, conceptualization, analysis, and writing.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work at the Florida International University was funded by the US Department of Energy, Basic Energy Sciences under the grant DE-FG02-04ER15570. AMM also acknowledges the Instructional & Research Computing Center (IRCC, web: https://ircc.fiu.edu) at the FIU for providing HPC computing resources that have contributed to the research results reported within this paper. Argonne is a U.S. Department of Energy laboratory managed by UChicago Argonne, LLC, under contract DE-AC02–06CH11357.

References

1. L. J. J. Catalan, B. Roberts and E. Rezaei, A low carbon methanol process using natural gas pyrolysis in a catalytic molten metal bubble reactor, Chem. Eng. J., 2023, 462, 142230.
2. A. Nadarajah, J. G. Lawrence and T. W. Hughes, Development and commercialization of vapor grown carbon nanofibers: A review, Key Eng. Mater., 2008, 380, 193–206.
3. R. Ciriminna, N. Zhang, M.-Q. Yang, F. Meneguzzo, Y.-J. Xu and M. Pagliaro, Commercialization of graphene-based technologies: a critical insight, Chem. Commun., 2015, 51, 7090–7095.
4. A. Dato, Graphene synthesized in atmospheric plasmas—A review, J. Mater. Res., 2019, 34, 214–230.
5. J. Diab, L. Fulcheri, V. Hessel, V. Rohani and M. Frenklach, Why turquoise hydrogen will Be a game changer for the energy transition, Int. J. Hydrogen Energy, 2022, 47, 25831–25848.
6. M. Faraday, The Chemical History of a Candle, Oxford University, New York, 2011.
7. M. Frenklach, A unifying picture of gas phase formation and growth of PAH, soot, diamond and graphite, in Carbon in the Galaxy: Studies From Earth and Space, ed. J. C. Tarter, S. Chang and D. J. DeFrees, NASA Conference Publication, 1990, 3061, pp. 259–273.
8. S. J. Klippenstein, From theoretical reaction dynamics to chemical modeling of combustion, Proc. Combust. Inst., 2017, 36, 77–111.
9. A. N. Morozov, A. M. Mebel and M. Frenklach, Acceleration of a chemical reaction due to nonequilibrium collisional dynamics: dimerization of polyaromatics, J. Phys. Chem. Lett., 2022, 13, 11528–11535.
10. M. Frenklach, Reaction mechanism of soot formation in flames, Phys. Chem. Chem. Phys., 2002, 4, 2028–2037.
11. M. Frenklach and A. M. Mebel, On the mechanism of soot nucleation, Phys. Chem. Chem. Phys., 2020, 22, 5314–5331.
12. J. W. Martin, M. Salamanca and M. Kraft, Soot inception: Carbonaceous nanoparticle formation in flames, Prog. Energy Combust. Sci., 2022, 88, 100956.
13. M. J. Thomson, Modeling soot formation in flames and reactors: recent progress and current challenges, Proc. Combust. Inst., 2023, 39, 805–823.
14. M. Frenklach and A. M. Mebel, Prenucleation chemistry of aromatics: a two-ring precursor?, Proc. Combust. Inst., 2023, 39, 825–833.
15. J. Ahmed and S. K. Mandal, Phenalenyl radical: Smallest polycyclic odd alternant hydrocarbon present in the graphene sheet, Chem. Rev., 2022, 122, 11369–11431.
16. K. Uchida and T. Kubo, Recent advances in the chemistry of phenalenyl, J. Synth. Org. Chem., 2016, 74, 1069–1077.
17. J. D. Bittner and J. B. Howard, Composition profiles and reaction mechanisms in a near-sooting premixed benzene/oxygen/argon flame, Symp. (Int.) Combust., 1981, 18, 1105–1116.
18. G. Vitiello, G. De Falco, F. Picca, M. Commodo, G. D'Errico, P. Minutolo and A. D'Anna, Role of radicals in carbon clustering and soot inception: a combined EPR and Raman spectroscopic study, Combust. Flame, 2019, 205, 286–294.
19. Z. D. Levey, B. A. Laws, S. P. Sundar, K. Nauta, S. H. Kable, G. da Silva, J. F. Stanton and T. W. Schmidt, PAH growth in flames and space: Formation of the phenalenyl radical, J. Phys. Chem. A, 2022, 126, 101–108.
20. L. Zhao, R. I. Kaiser, W. Lu, M. Ahmed, A. D. Oleinikov, V. N. Azyazov, A. M. Mebel, A. H. Howlader and S. F. Wnuk, Gas phase formation of phenalene via 10π-aromatic, resonantly stabilized free radical intermediates, Phys. Chem. Chem. Phys., 2020, 22, 15381–15388.
21. D. P. Porfiriev, V. N. Azyazov and A. M. Mebel, Conversion of acenaphthalene to phenalene via methylation: a theoretical study, Combust. Flame, 2020, 213, 302–313.
22. Z. Mou, K. Uchida, T. Kubo and M. Kertesz, Evidence of σ- and π-dimerization in a series of phenalenyls, J. Am. Chem. Soc., 2014, 136, 18009–18022.
23. L.-A. Lieske, M. Commodo, J. W. Martin, K. Kaiser, V. Benekou, P. Minutolo, A. D’Anna and L. Gross, Portraits of soot molecules reveal pathways to large aromatics, five-/seven-membered rings, and inception through π-radical localization, ACS Nano, 2023, 17, 13563–13574.
24. A. Schneiker, I. P. Csonka and G. Tarczay, Hydrogenation and dehydrogenation reactions of the phenalenyl radical/1H-phenalene system at low temperatures, Chem. Phys. Lett., 2020, 743, 137183.
25. K.-H. Homann, Fullerenes and soot formation—New pathways to large particles in flames, Angew. Chem., Int. Ed., 1998, 37, 2434–2451.
26. A. Keller, R. Kovacs and K. H. Homann, Large molecules, ions, radicals and small soot particles in fuel-rich hydrocarbon flames. Part IV. Large polycyclic aromatic hydrocarbons and their radicals in a fuel-rich benzene–oxygen flame, Phys. Chem. Chem. Phys., 2000, 2, 1667–1675.
27. A. D. Becke, Density-functional thermochemistry. III. The role of exact exchange, J. Chem. Phys., 1993, 98, 5648–5652.
28. C. T. Lee, W. T. Yang and R. G. Parr, Development of the Colle-Salvetti correlation-energy formula into a functional of the electron density, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785–789.
29. A. G. Baboul, L. A. Curtiss, P. C. Redfern and K. Raghavachari, Gaussian-3 theory using density functional geometries and zero point energies, J. Chem. Phys., 1999, 110, 7650–7657.
30. L. A. Curtiss, K. Raghavachari, P. C. Redfern, A. G. Baboul and J. A. Pople, Gaussian-3 theory using coupled cluster energies, Chem. Phys. Lett., 1999, 314, 101–107.
31. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson and H. Nakatsuji, et al., Gaussian 16 (Revision C.01), Gaussian Inc., Wallingford, CT, 2019.
32. H.-J. Werner, P. J. Knowles, P. Celani, W. Györffy, A. Hesselmann, D. Kats, G. Kinizia, A. Köhn, T. Korona and D. Kreplin, et al., MOLPRO (2021.2), a package of ab initio programs, see https://www.molpro.net, 2021.
33. H. J. Werner and P. J. Knowles, A second order multiconfiguration SCF procedure with optimum convergence, J. Chem. Phys., 1985, 82, 5053–5063.
34. P. J. Knowles and H.-J. Werner, An efficient second-order MC SCF method for long configuration expansions, Chem. Phys. Lett., 1985, 115, 259–267.
35. J.-D. Chai and M. Head-Gordon, Long-range corrected hybrid density functionals with damped atom-atom dispersion corrections, Phys. Chem. Chem. Phys., 2008, 10, 6615–6620.
36. A. W. Jasper, NST, A Spin-Forbidden Nonadiabatic Flux Code, https://tcg.cse.anl.gov/papr/codes/nst.html, (accessed August 30, 2023).
37. J. N. Harvey and M. Aschi, Spin-forbidden dehydrogenation of methoxy cation: a statistical view, Phys. Chem. Chem. Phys., 1999, 1, 5555–5563.
38. J. N. Harvey, Understanding the kinetics of spin-forbidden chemical reactions, Phys. Chem. Chem. Phys., 2007, 9, 331–343.
39. A. W. Jasper, Multidimensional effects in nonadiabatic statistical theories of spin-forbidden kinetics: a case study of 3O + CO → CO₂, J. Phys. Chem. A, 2015, 119, 7339–7351.
40. Y. Georgievskii, J. A. Miller, M. P. Burke and S. J. Klippenstein, Reformulation and solution of the master equation for multiple-well chemical reactions, J. Phys. Chem. A, 2013, 117, 12146–12154.
41. Y. Georgievskii and S. J. Klippenstein, Master Equation System Solver (MESS) (2016.3.23), 2015.
42. J. Troe, Theory of thermal unimolecular reactions at low pressures. I. Solutions of the master equation, J. Chem. Phys., 1977, 66, 4745–4757.
43. A. W. Jasper, C. M. Oana and J. A. Miller, “Third-Body” collision efficiencies for combustion modeling: Hydrocarbons in atomic and diatomic baths, Proc. Combust. Inst., 2015, 35, 197–204.
44. A. M. Mebel, Y. Georgievskii, A. W. Jasper and S. J. Klippenstein, Temperature- and pressure-dependent rate coefficients for the HACA pathways from benzene to naphthalene, Proc. Combust. Inst., 2017, 36, 919–926.
45. A. M. Mebel, Y. Georgievskii, A. W. Jasper and S. J. Klippenstein, Pressure-dependent rate constants for PAH growth: formation of indene and its conversion to naphthalene, Faraday Discuss., 2016, 195, 637–670.
46. M. Frenklach, R. I. Singh and A. M. Mebel, On the low-temperature limit of HACA, Proc. Combust. Inst., 2019, 37, 969–976.
47. T. Yang, R. I. Kaiser, T. P. Troy, B. Xu, O. Kostko, M. Ahmed, A. M. Mebel, M. V. Zagidullin and V. N. Azyazov, HACA's heritage: a free-radical pathway to phenanthrene in circumstellar envelopes of asymptotic giant branch stars, Angew. Chem., Int. Ed., 2017, 56, 4515–4519.
48. L. B. Tuli and A. M. Mebel, Formation of phenanthrene via H-assisted isomerization of 2-ethynylbiphenyl produced in the reaction of phenyl with phenylacetylene, Int. J. Chem. Kinet., 2020, 52, 875–883.
49. A. S. Savchenkova, I. V. Chechet, S. G. Matveev, M. Frenklach and A. M. Mebel, Formation of phenanthrenyl radicals via the reaction of acenaphthyl with acetylene, Proc. Combust. Inst., 2021, 38, 1441–1448.
50. G. R. Galimova, I. A. Medvedkov and A. M. Mebel, The role of methylaryl radicals in the growth of polycyclic aromatic hydrocarbons: the formation of five-membered rings, J. Phys. Chem. A, 2022, 126, 1233–1244.
51. A. D. Abid, J. Camacho, D. A. Sheen and H. Wang, Quantitative measurement of soot particle size distribution in premixed flames—The burner-stabilized stagnation flame approach, Combust. Flame, 2009, 156, 1862–1870.
52. F. Xu, P. B. Sunderland and G. M. Faeth, Soot formation in laminar premixed ethylene/air flames at atmospheric pressure, Combust. Flame, 1997, 108, 471–493.
53. G. P. Smith, Y. Tao and H. Wang, Foundational fuel chemistry model Version 1.0 (FFCM-1), 2016, https://web.stanford.edu/group/haiwanglab/FFCM1/pages/FFCM1.html, (accessed August 18, 2023).
54. D. G. Goodwin, H. K. Moffat, I. Schoegl, R. L. Speth and B. W. Weber, Cantera: An object-oriented software toolkit for chemical kinetics, thermodynamics, and transport processes (Version 2.6.0), 2022 DOI:10.5281/zenodo.6387882.
55. Matlab, https://www.mathworks.com/, (accessed September 2, 2023).
56. J. Appel, H. Bockhorn and M. Frenklach, Kinetic modeling of soot formation with detailed chemistry and physics: laminar premixed flames of C₂ hydrocarbons, Combust. Flame, 2000, 121, 122–136.
57. A. M. Valencia-López, F. Bustamante, A. Loukou, B. Stelzner, D. Trimis, M. Frenklach and N. A. Slavinskaya, Effect of benzene doping on soot precursors formation in non-premixed flames of producer gas (PG), Combust. Flame, 2019, 207, 265–280.
58. D. T. Gillespie, Exact stochastic simulation of coupled chemical reactions, J. Phys. Chem., 1977, 81, 2340–2361.
59. R. Whitesides and M. Frenklach, Detailed kinetic Monte Carlo simulations of graphene-edge growth, J. Phys. Chem. A, 2010, 114, 689–703.
60. Y. Carissan, D. Hagebaum-Reignier, N. Prcovic, C. Terrioux and A. Varet, How constraint programming can help chemists to generate benzenoid structures and assess the local aromaticity of benzenoids, Constraints, 2022, 27, 192–248.
61. L. B. Tuli, A. M. Mebel and M. Frenklach, Bay capping via acetylene addition to polycyclic aromatic hydrocarbons: mechanism and kinetics, Proc. Combust. Inst., 2023, 39, 969–977.
62. V. V. Kislov, N. I. Islamova, A. M. Kolker, S. H. Lin and A. M. Mebel, Hydrogen Abstraction Acetylene Addition and Diels–Alder mechanisms of PAH formation:  a detailed study using first principles calculations, J. Chem. Theor. Comput., 2005, 1, 908–924.
63. A. S. Semenikhin, A. S. Savchenkova, I. V. Chechet, S. G. Matveev, Z. Liu, M. Frenklach and A. M. Mebel, Rate constants for H abstraction from benzo(a)pyrene and chrysene: a theoretical study, Phys. Chem. Chem. Phys., 2017, 19, 25401–25413.
64. B. D. Adamson, S. A. Skeen, M. Ahmed and N. Hansen, Detection of aliphatically bridged multi-core polycyclic aromatic hydrocarbons in sooting flames with atmospheric-sampling high-resolution tandem mass spectrometry, J. Phys. Chem. A, 2018, 122, 9338–9349.
65. B. A. Adamson, S. A. Skeen, M. Ahmed and N. Hansen, Nucleation of soot: experimental assessment of the role of polycyclic aromatic hydrocarbon (PAH) dimers, Z. Phys. Chem., 2020, 234, 1295–1310.
66. H. Jin, J. Hao, J. Yang, J. Guo, Y. Zhang, C. Cao and A. Farooq, Experimental and kinetic modeling study of α-methyl-naphthalene pyrolysis: Part II. PAH formation, Combus. Flame, 2021, 233, 111530.
67. K. O. Johansson, M. P. Head-Gordon, P. E. Schrader, K. R. Wilson and H. A. Michelsen, Resonance-stabilized hydrocarbon-radical chain reactions may explain soot inception and growth, Science, 2018, 361, 997–1000.
68. M. Frenklach, D. W. Clary, W. C. Gardiner, Jr. and S. E. Stein, Detailed kinetic modeling of soot formation in shock-tube pyrolysis of acetylene, Proc. Combust. Inst., 1985, 20, 887–901.
69. M. Weissman and S. W. Benson, Pyrolysis of methyl chloride, a pathway in the chlorine-catalyzed polymerization of methane, Int. J. Chem. Kinet., 1984, 16, 307–333.
70. E. Georganta, R. K. Rahman, A. Raj and S. Sinha, Growth of polycyclic aromatic hydrocarbons (PAHs) by methyl radicals: pyrene formation from phenanthrene, Combust. Flame, 2017, 185, 129–141.
71. B. Shukla, A. Miyoshi and M. Koshi, Role of methyl radicals in the growth of PAHs, J. Am. Soc. Mass Spectrom., 2010, 21, 534–544.
72. N. Hansen, M. Schenk, K. Moshammer and K. Kohse-Höinghaus, Investigating repetitive reaction pathways for the formation of polycyclic aromatic hydrocarbons in combustion processes, Combust. Flame, 2017, 180, 250–261.
73. V. S. Krasnoukhov, D. P. Porfiriev, I. P. Zavershinskiy, V. N. Azyazov and A. M. Mebel, Kinetics of the CH₃ + C₅H₅ reaction: a theoretical study, J. Phys. Chem. A, 2017, 121, 9191–9200.
74. A. M. Mebel and M. Frenklach, Cleavage of an aromatic ring and radical migration, Faraday Discuss., 2022, 238, 512–528.
75. M. Frenklach, Z. Liu, R. I. Singh, G. R. Galimova, V. N. Azyazov and A. M. Mebel, Detailed, sterically-resolved modeling of soot oxidation: role of O atoms, interplay with particle nanostructure, and emergence of inner particle burning, Combust. Flame, 2018, 188, 284–306.
76. M. Frenklach, Detailed modeling of soot particles nucleation and growth. Poster, presented at the Twenty-Second International Symposium on Combustion, Seattle, Washington, August 14-19, 1988.
77. M. Frenklach, Comment, in Soot Formation in Combustion: An International Round Table Discussion, ed. H. Jander and H. G. Wagner, Göttingen, 1989, pp. 152–153.
78. M. Frenklach and H. Wang, Detailed modeling of soot particle nucleation and growth, Proc. Combust. Inst., 1991, 23, 1559–1566.
79. M. Frenklach and W. C. Gardiner, Representation of multistage mechanisms in detailed computer modeling of polymerization kinetics, J. Phys. Chem., 1984, 88, 6263–6266.
80. M. Frenklach, Computer modeling of infinite reaction sequences: a chemical lumping, Chem. Eng. Sci., 1985, 40, 1843–1849.
81. M. Frenklach and H. Wang, Detailed mechanism and modeling of soot particle formation, in Soot Formation in Combustion: Mechanisms and Models, ed. H. Bockhorn, Springer-Verlag, Heidelberg, 1994, pp. 165–192.
82. M. Frenklach and S. J. Harris, Aerosol dynamics modeling using the method of moments, J. Colloid Interface Sci., 1987, 118, 252–261.
83. J. H. Miller, K. C. Smyth and W. G. Mallard, Calculations of the dimerization of aromatic hydrocarbons: implications for soot formation, Proc. Combust. Inst., 1985, 20, 1139–1147.
84. H. Sabbah, L. Biennier, S. J. Klippenstein, I. R. Sims and B. R. Rowe, Exploring the role of PAHs in the formation of soot: pyrene dimerization, J. Phys. Chem. Lett., 2010, 1, 2962–2967.
85. P. Elvati and A. Violi, Thermodynamics of poly-aromatic hydrocarbon clustering and the effects of substituted aliphatic chains, Proc. Combust. Inst., 2013, 34, 1837–1843.
86. T. S. Totton, A. J. Misquitta and M. Kraft, A quantitative study of the clustering of polycyclic aromatic hydrocarbons at high temperatures, Phys. Chem. Chem. Phys., 2012, 14, 4081–4094.
87. A. S. Semenikhin, A. S. Savchenkova, I. V. Chechet, S. G. Matveev, M. Frenklach and A. M. Mebel, On the mechanism of soot nucleation. II. E-bridge formation at the PAH bay, Phys. Chem. Chem. Phys., 2020, 22, 17196–17204.
88. M. Frenklach, A. S. Semenikhin and A. M. Mebel, On the mechanism of soot nucleation. III. The fate and facility of the E-bridge, J. Phys. Chem. A, 2021, 125, 6789–6795.
89. H. Wang, Formation of nascent soot and other condensed-phase materials in flames, Proc. Combust. Inst., 2011, 33, 41–67.
90. P. K. Selvakumar, J. W. Martin, M. D. Lorenzo, M. Paskevicius and C. E. Buckley, Role of π-radical localization on thermally stable cross-links between polycyclic aromatic hydrocarbons, J. Phys. Chem. A, 2023, 127, 6945–6952.
91. L. Pascazio, J. W. Martin, A. Menon, D. Hou, X. You and M. Kraft, Aromatic penta-linked hydrocarbons in soot nanoparticle formation, Proc. Combust. Inst., 2021, 38, 1525–1532.
92. J. W. Martin, L. Pascazio, A. Menon, J. Akroyd, K. Kaiser, F. Schulz, M. Commodo, A. D’Anna, L. Gross and M. Kraft, π-diradical aromatic soot precursors in flames, J. Am. Chem. Soc., 2021, 143, 12212–12219.
93. A. Menon, J. Martin, G. Leon, D. Hou, L. Pascazio, X. You and M. Kraft, Reactive localized π-radicals on rim-based pentagonal rings: properties and concentration in flames, Proc. Combust. Inst., 2021, 38, 565–573.
94. J. Huang and M. Kertesz, Theoretical analysis of intermolecular covalent π–π bonding and magnetic properties of phenalenyl and spiro-biphenalenyl radical π-dimers, J. Phys. Chem. A, 2007, 111, 6304–6315.
95. A. Das, T. Müller, F. Plasser and H. Lischka, Polyradical character of triangular non-Kekulé structures, zethrenes, p-quinodimethane-linked bisphenalenyl, and the Clar goblet in comparison: An extended multireference study, J. Phys. Chem. A, 2016, 120, 1625–1636.
96. M. Frenklach and A. M. Mebel, On the mechanism of soot nucleation. IV. Molecular growth of the flattened E-bridge, J. Phys. Chem. A, 2022, 126, 9259–9267.
97. M. L. Botero, E. M. Adkins, S. González-Calera, H. Miller and M. Kraft, PAH structure analysis of soot in a non-premixed flame using high-resolution transmission electron microscopy and optical band gap analysis, Combust. Flame, 2016, 164, 250–258.
98. M. L. Botero, Y. Sheng, J. Akroyd, J. Martin, J. A. H. Dreyer, W. Yang and M. Kraft, Internal structure of soot particles in a diffusion flame, Carbon, 2019, 141, 635–642.
99. R. S. Jacobson, A. R. Korte, A. Vertes and J. H. Miller, The molecular composition of soot, Angew. Chem., Int. Ed., 2020, 59, 4484–4490.
100. H. Wang and M. Frenklach, Transport properties of polycyclic aromatic hydrocarbons for flame modeling, Combust. Flame, 1994, 96, 163–170.
101. S. H. Bauer, Comment, in Particulate Carbon: Formation During Combustion, ed. D. C. Siegla and G. W. Smith, Plenum, New York, 1981, pp. 169–173.
102. M. Frenklach, Kinetics of powder nucleation, in Ceramic Powder Science III, ed. G. L. Messing, S.-I. Hirano and H. Hausner, American Ceramic Society, Westerville, Ohio, 1990, pp. 239–250.

---

Footnote

† Electronic supplementary information (ESI) available: Reaction rate constants and MESS input files. See DOI: https://doi.org/10.1039/d4cp00096j

---