Reaction pathways for growth of polycyclic aromatic hydrocarbons under combustion conditions, a DFT study

Models of soot formation include the growth of polycyclic aromatic hydrocarbons (PAHs) as an essential step.

We investigate two reaction pathways for growth of PAHs, both initiated by radical formation.

The first type, the so called HACA mechanism, is based on stepwise acetylene addition and we identify promising reactions not considered before.

Examples are the synthesis of perylene and phenanthro[1,10,9,8–opqra]perylene.

The second type of reaction concerns direct condensation of phenyl, naphthyl or anthracenyl with another PAH.

Pathways for these reactions require merely up to 100 kJ mol−1 for the activation barrier, which makes them viable candidates for PAH growth.


Soot is a typical byproduct of hydrocarbon combustion; this is undesirable since it indicates incomplete combustion and also entails a considerable health hazard.1

Soot is at the same time a viable material if prepared with well defined properties, that can be achieved in pyrolysis from, e.g. methane, acetylene or other hydrocarbons.

So produced, it has many applications in chemical industry and material science.2–5

Detailed understanding of mechanisms of soot formation in combustions, pyrolysis, and also techniques such as carbon vapor deposition (CVD) are essential for controlling these processes.

It has been the main aim of the present work to better characterize the elementary reactions in the gas phase, which play a role in soot formation.

The central issue is here the formation of polycyclic aromatic hydrocarbons (PAHs), which are the main constituents of soot; the PAHs are also health hazards since they are carcinogenic.1

The formation of the first aromatic ring under pyrolysis conditions has been studied by various groups.

It appears that there are different types of growth reactions pathways depending on the fuel and actual conditions.6–9

Starting from methane or acetylene as precursors Frenklach et al. propose a reaction for the formation of the first aromatic ring from n-C4H3 and acetylene.10,11

Bittner and Howard suggest in addition the reaction of n-C4H5 with acetylene occurring especially at low temperatures.

Miller and Melius pointed to a problem in these C4-routes: n-C4H3 and n-C4H5 should convert to the corresponding resonance stabilized isomers (iso-C4H3 and iso-C4H5).12

Therefore, they recommend together with others6,13,14 a C3-route over propargyl radicals.

For further growth (beyond benzene), the most popular mechanism is the hydrogen-abstraction/acetylene-addition, which is explicitly suggested by Bockhorn et al.15

Independently Frenklach and coworkers refined a similar sequence of reactions10 and termed it the ‘HACA’-mechanism (Hydrogen-Abstraction-C2H2-Addition).16

Böhm et al. have proposed a combination of a ring–ring condensation and acetylene addition,17 since reactions with acetylene alone apparently cannot describe the timescale of soot formation.

A related model is proposed by Dong and Hüttinger, the so called particle-filler model:18 aromatic hydrocarbons constitute the molecular particles and acetylene the molecular filler.

In this model the initial step leads to a biaryl bond, then followed by a dehydrocyclization reaction forming a five- or six-membered ring.

Bays of PAHs involving three or four carbon atoms can be closed by addition of acetylene leading to five- and six-membered rings.

With their model18 they are able to predict the texture of pyrolytic carbon: at low temperatures, low pressures, short residence times and high surface area/volume ratios, formation of sufficient amount of PAHs is limited thus leading to low-ordered carbon.

On the other hand one gets a highly ordered carbon under high temperatures, high pressures, long residence times and low surface area/volume ratios.

However, the temperatures should not be too high to avoid too fast growth.18

The present computational study is carried out to obtain energetic data for some reactions of interest in soot formation, as indicated above.

We have shown previously that semiempirical PM3 calculations give occasionally considerable errors, e.g. 100 kJ mol−1 for the addition of a hydrogen atom to pyrene.19

In this work we describe different reaction pathways leading to PAHs; first we determine barriers and reaction energies for the HACA mechanism for naphthalene, acenaphthylene and phenanthrene as examples.

As far as PM3 data are available,20 we compare these with our calculations.

We follow here the HACA-mechanism and do not consider reactions like PAH + C2H; since the latter require additional hydrogen migration and can only proceed by a perturbation of the aromatic system.

Second we study aromatic condensation of PAHs; this should have a bearing on the model of Dong and Hüttinger.18

We assume that, as in the HACA process, condensation is initiated by radical formation, and then consider reactions of phenyl, naphthyl and anthracenyl with another PAH to yield for example fluoranthene or perylene.

Reaction pathways of the polyaromatic growth have clearly been considered before, but these investigations have either been carried out at the PM3 level,20 or are limited to systems with ten carbon atoms at most;21–23 aromatic condensations have not been modeled before to our knowledge.

Computational details

Geometry optimizations have been carried out to locate minima and transition states with gradient corrected DFT methods as implemented in the TURBOMOLE program package.24–26

We use mainly the functional B3LYP,27,28 since this is the standard tool for these purposes.

In addition, all investigations have been done using the BP86 functional29–32 in combination with the very efficient RI-J approximation,33,34 which reduces the CPU times to about 10% of corresponding B3LYP calculations.

This was done to ascertain the reliability of BP86, which offers the prospect to treat much larger systems.

The standard TURBOMOLE basis set SV(P) (split valence augmented by a polarization d set on C) is employed35,36 and test calculations are performed with larger basis sets like TZVPP,37 and restricted open shell (ROHF-) RIMP2,38,39 and also ROHF-RICC.240

Redundant internal coordinates are used in the structure optimization of minima.41

Optimization of saddle points is performed with the TRIM (trust radius image minimization)42 method now implemented in TURBOMOLE.

Harmonic frequency calculations are carried out for all species;43,44 all minima presented show only real frequencies and the first-order saddle points have only a single imaginary frequency.

Reaction enthalpies are computed within the harmonic approximation using frequencies obtained from analytical second derivatives at B3LYP/SV(P) and BP86/SV(P) level, respectively.

Accuracy of the method

Let us first point out that accuracy requirements are relatively moderate for present purposes, since we are dealing with high-temperature reactions; pyrolysis is typically carried out at 1100 K–1300 K. To estimate the error in the barrier energies we consider the hydrogen abstraction with a hydrogen atom from the gas phase, since this is a crucial step for many of the reactions considered below.

As first example we chose a radical reaction with small molecules for which experimental data are available CH4 + H → CH3HH → CH3 + H2.The barrier of the transition state calculated on DFT/B3LYP/SV(P) is +38 kJ mol−1 (+25 kJ mol−1 for DFT/BP86/SV(P)) and the experimental result is +50 kJ mol−1.45

The next system considered is the reaction C6H6 + H → C6H5 + H2.

In Table 1 we compare various methods and basis sets.

ΔE1 denotes the difference between educts and transition state (barriers) and ΔE2 the difference between educts and products (reaction energies).

Here, one can see that the computed ΔE decrease with increasing basis, and the smallest basis set gives values near experimental ΔE2, particularly for DFT/B3LYP/SV(P); this has also been pointed out by Houk et al.46

The large deviation of the ROHF-RIMP2 and ROHF-RICC2 methods is caused by the (partial) multireference character of phenyl;47 the calculated values are simply not reliable.

We have further considered the next larger system, naphthalene + hydrogen, the results are very similar to those listed in Table 1.

B3LYP typically underestimates barriers by 10–35%,48,49 which is a well known and presently unavoidable deficiency, at least for systems of this size.

We follow previous investigations20,50 and report ΔHR values under standard conditions within the harmonic approximation, we will comment on this in the summary.

Results and discussion

Acetylene addition

Frenklach and Wang have proposed that the growth of soot particles on surfaces takes place via the so called HACA mechanism.16

This consists of the following steps: hydrogen abstraction creating a surface radical, acetylene addition to the radical formed, and ring cyclization as shown in Fig. 1.

The first reaction considered in this work, subsection 1, is the generation of acenaphthylene [6] from naphthyl radical [1] and acetylene, a five-ring cyclization at the smallest system representing a zigzag border.47

The second reaction, subsection 2, is a transformation of the five-membered ring to a six-membered ring by adding acetylene to an acenaphthyl radical [6a] resulting in a phenanthryl radical [17].

This reaction was proposed by Frenklach et al20. and we will compare their PM3 results with DFT (BP86 and B3LYP).

The next reaction, subsection 3, will also be compared with PM3: building a six-membered ring on phenanthryl [17] to form pyrene [23].

In subsection 4 we compare the HACA mechanism of Frenklach and the Bittner–Howard mechanism,51 where two acetylene are added in sequence to acenaphthyl [6a] leading to aceanthrylene [31].

C10H7 [1] + C2H2 → C12H8 [6] + H

Fig. 2. shows the reaction enthalpies at 298 K for the acetylene addition to naphthyl radical [1] to yield acenaphthylene [6] and hydrogen.

In Table 2 we compare calculated B3LYP and BP86 results for elementary steps of this reaction.

The reaction is exothermic by 163 kJ mol−1 for B3LYP (181 kJ mol−1 for BP86), the highest barrier of the total reaction path (relative to the educts) is 19 kJ mol−1 for B3LYP (8 kJ mol−1 for BP86).

The reaction begins with the formation of [2] which proceeds to either [3] via hydrogen migration or to [4] under hydrogen atom loss.

At B3LYP level the way to [3] shows a barrier of 26 kJ mol−1 (4 kJ mol−1 BP86), significantly lower than the one to [4] of 150 kJ mol−1 (147 kJ mol−1 BP86).

However, the reaction channel to [4] is exothermic and the barrier is moderate considering typical reaction conditions of pyrolysis.52

The barrier of the ring cyclization of [3] to [5] is 42 kJ mol−1 for B3LYP (28 kJ mol−1 for BP86).

Forming [6] under hydrogen loss has an activation barrier of 192 kJ mol−1, the transition state is just 2 kJ mol−1 above the energy of the final product (BP86 shows no additional barrier).

C12H7 [6a] + C2H2 → C14H9 [17]

The second reaction, Fig. 3 and Table 3, is the transformation from a five- to a six-membered ring.

The reaction is exothermic by 312 kJ mol−1 for B3LYP (320 kJ mol−1 for BP86).

The HACA mechanism proceeds through formation of [11] and [16], respectively and the highest barrier of the reaction path is found to be 227 kJ mol−1 for B3LYP (227 kJ mol−1 for BP86) for the formation of [11].

The treatment of the first steps of the reaction path up to [10] is straightforward.

But then we run into a problem: Educt [6a] and product [17] are radicals with clear cut doublet character (one singly occupied molecular orbital); a visualization of spin densities further shows the unpaired spin to be localized on a single center.

The intermediates [11] and [16] give rise to three radical centers, which have to be coupled to a doublet.

This is not possible in a clear cut way with standard DFT.

Since the three singly occupied molecular orbitals in [11] are located on different centers, as indicated in Fig. 3, one has small exchange integrals and energies obtained with different spin distributions (ααα vs. ααβ etc.) are virtually identical.

The differences amount to merely 4 kJ mol−1 for B3LYP (and also for BP86).

(On the DFT level [16] is actually a minimum for the quartet state, the energy in Table 3 is given in parentheses).

Our investigation leads to a new pathway starting from [10] over TS12, [14] and TS14 to yield [15], solid black line in Fig. 3.

This proceeds over clear cut doublets only (〈S2〉 = 0.78 for TS12) and avoids the problems just discussed.

The barrier from [10] to TS12 is 182 kJ mol−1 for B3LYP (153 kJ mol−1 BP86).

This should be the favored reaction pathway.

A comparison with Frenklach’s PM3 results20 shows reasonable agreement with present results up to [10], but we do not find indications that intermediates with three radical centers play a role for this reaction.

C14H9 [17] + C2H2 → C16H10 [23] + H

The six-membered ring cyclization is depicted in Fig. 4 and the energetic results are presented in Table 4.

Here also PM3 results are available.20

We find an exothermic reaction of 257 kJ mol−1 for B3LYP (268 kJ mol−1 for BP86).

The highest barrier of the total path is calculated to be 29 kJ mol−1 for B3LYP (15 kJ mol−1 for BP86).

Starting after the formation of [18] the radical has three possibilities to react.

The first way is described in the two reactions discussed above: under hydrogen loss we get [21], which presents the same behavior as described in subsection 3.1.2 (C12H7 [6a] + C2H2 → C14H9 [17]).

PM3 results show an endothermic reaction and DFT an exothermic pathway.

The second way leads to [19] via ring cyclization and the last way is a H-migration from the C2H2 group to the six-membered ring [20] followed by ring cyclization [22].

Both radicals, [19] and [22], decompose to pyrene [23].

Discrepancies between PM3 and B3LYP (or BP86) are considerable, Table 4: the barrier height of TS20 from PM3 is 100 kJ mol−1, while it is 16 kJ mol−1 for DFT/B3LYP (14 kJ mol−1 for BP86).

Furthermore the reaction enthalpies for the formation of pyrene differ by about 107 kJ mol−1 between B3LYP and PM3 (118 kJ mol−1 between BP86 and PM3).

This has also been pointed out in a previous paper19.

C12H7 [6a] + 2 C2H2 → C16H10 [31] + H

For this reaction we compare Frenklach’s HACA mechanism (solid line in Fig. 5), and the mechanism of Bittner and Howard,51 (dashed line in Fig. 5).

Bauschlicher and coworkers have studied both mechanisms for a similar reaction: benzene + 2 C2H2 → naphthalene.53

They found that both mechanisms have low barriers.

The reaction is exothermic with 449 kJ mol−1 for B3LYP (477 kJ mol−1 for BP86, see Table 5).

The first step is the exothermic formation of [24] from [6a] and acetylene, which shows only a small barrier of 18 kJ mol−1 for B3LYP (7 kJ mol−1 for BP86) and this is the highest barrier of the total reaction pathway.

From this point the two pathways differ because in the HACA mechanism the second acetylene is added later than in the Bittner–Howard mechanism.

We first consider the Frenklach mechanism (solid line in Fig. 5).

The hydrogen migration from the six-membered ring to the C2H2 group [26] has a barrier of 104 kJ mol−1 for B3LYP (84 kJ mol−1 for BP86).

Again an acetylene addition leads to [28] with a barrier of 25 kJ mol−1 for B3LYP (11 kJ mol−1 for BP86).

A ring cyclization [30] takes place without any barrier.

The reaction ends forming aceanthrylene [31] by decomposition of [30].

The Bittner–Howard mechanism on [24] is a direct addition of acetylene at the C2H2 group of [25] with a barrier of 20 kJ mol−1 for B3LYP (13 kJ mol−1 for BP86).

This barrier is significantly lower than the step from [24] to [26].

The 1,6-H-shift [27] has no barrier for BP86 calculations and a small one with the B3LYP functional.

Another small activation barrier from [27] to [29] is found for the ring cyclization.

Decomposition of [29] leads to aceanthrylene [31].

As a summary, we find that the Bittner–Howard mechanism is slightly favored, since the barrier for [24] to [25] is only 23 kJ mol−1 for B3LYP (13 kJ mol−1 for BP86).

This result is in line with that of Bauschlicher et al53. for the analogous reaction of benzene and two acetylene to form naphthalene.

The formation of [9] under hydrogen loss is not further investigated because the other pathways are significantly lower.

Aromatic condensation

In the second part we consider the condensation of aromatic compounds.

Two PAHs are connected to form a larger PAH under the formation of hydrogen molecule(s).

We start again with number [1] counting the molecules.

In the first reaction, subsection 1, benzene [1] reacts with naphthyl radical [2a] forming fluoranthene [8].

The next reaction (subsection 2), is the condensation of naphthalene [2] with naphthyl radical [2a] to perylene [13].

In the last process, subsection 3, anthracene [14] reacts with an anthracenyl radical [14a] to form phenanthro[1,10,9,8–opqra]perylene (phenanthroperylene) [22].

C6H6 [1] + C10H7 [2a] → C16H10 [8] + H2 + H

Fig. 6 and Table 6 show the aromatic condensation of benzene [1] with a naphthyl radical [2a].

The reaction is endothermic by 24 kJ mol−1 for B3LYP (14 kJ mol−1 for BP86) and the highest barrier of the reaction path, relative to the educts, is obtained to be 85 kJ mol−1 for B3LYP (60 kJ mol−1 for BP86).

The analogous reaction of phenyl [1a] and naphthalene [2] is also calculated (dashed line).

[3] is lower in energy than [4] by 39 kJ mol−1 for B3LYP (37 kJ mol−1 for BP86).

This is due to the stabilization of the radical through the π-system.

The next step is the separation of a hydrogen atom leading to [5].

The energy of the barrier from [3] to [5] is found to be 132 kJ mol−1 for B3LYP (126 kJ mol−1 for BP86) and from [4] to [5] it is 104 kJ mol−1 for B3LYP (96 kJ mol−1 BP86).

The following step is an assisted hydrogen abstraction under formation of H2 and [6] with a barrier of 51 kJ mol−1 for B3LYP (36 kJ mol−1 for BP86).

The barrier of the ring cyclization forming a fluoranthenyl radical [7] has been computed to be 52 kJ mol−1 for B3LYP (37 kJ mol−1 for BP86).

The final step is the decomposition to fluoranthene [8], through a barrier of 113 kJ mol−1 for B3LYP (108 kJ mol−1 for BP86).

C10H7 [2a] + C10H8 [2] → C20H12 [13] + H2 + H

The second example for the aromatic condensation is the reaction of naphthalene [2] with a naphthyl radical [2a] to perylene [13], see Fig. 7.

This reaction is slightly exothermic by 5 kJ mol−1 for B3LYP (19 kJ mol−1 for BP86) and the highest barrier of the total reaction pathway is calculated to be 65 kJ mol−1 for B3LYP (45 kJ mol−1 for BP86).

The first barrier, to get [9], is 27 kJ mol−1 for B3LYP (16 kJ mol−1 for BP86) (Table 7).

The radical is again stabilized through the π-system, as mentioned above for [3].

Hydrogen loss leads to [10] over a barrier of 131 kJ mol−1 for B3LYP (126 kJ mol−1 for BP86).

The following step is a hydrogen abstraction by forming a hydrogen molecule and [11].

The barrier is 68 kJ mol−1 for B3LYP (50 kJ mol−1 for BP86).

The ring cyclization [12] has a small barrier of 25 kJ mol−1 for B3LYP (13 kJ mol−1 for BP86).

Separation of a hydrogen atom gives perylene [13].

The barrier is 104 kJ mol−1 for B3LYP (108 kJ mol−1 for BP86).

C14H9 [14a] + C14H10 [14] → C28H14 [22] + 2 H2 + H

The reaction of anthracene [14] with an anthracenyl radical [14a] is the third aromatic condensation, see Fig. 8 and Table 8.

The reaction is exothermic by 24 kJ mol−1 for B3LYP (47 kJ mol−1 for BP86) and the highest barrier of the total reaction pathway, relative to the educts, is calculated to be 97 kJ mol−1 for B3LYP (62 kJ mol−1 for BP86).

Starting with the addition of the two molecules one gets [15].

The computed barrier is 21 kJ mol−1 for B3LYP (14 kJ mol−1 for BP86).

Hydrogen separation leading to [16] has a barrier of 147 kJ mol−1 for B3LYP (139 kJ mol−1 for BP86).

The following H abstraction under formation of H2 and [17] has a barrier height of 59 kJ mol−1 for B3LYP (42 kJ mol−1 for BP86) similar to the two aromatic condensations discussed in subsection 1 and 2.

A first ring cyclization follows, introducing a steric strain in the molecule [18].

The barrier is 17 kJ mol−1 for B3LYP (6 kJ mol−1 for BP86).

The next step is the second loss of hydrogen to [19].

Then follows another hydrogen abstraction with a similar barrier as for the first one: 51 kJ mol−1 for B3LYP (26 kJ mol−1 for BP86) for [20] and H2.

The second ring cyclization leading to [21] occurs next.

At BP86 level no barrier is obtained; B3LYP leads to an insignificant barrier of 1 kJ mol−1.

Consequently, the ring closure takes place directly.

The last step is a hydrogen loss with a barrier of 106 kJ mol−1 for B3LYP (102 kJ mol−1 for BP86).

Brief discussion

We would like to mention briefly common features and differences of the aromatic condensation reactions.

The crucial step for the smaller cases, see Fig. 6 and Table 6, is clearly the ring closure, e.g. from [6] to [7].

The corresponding barrier is reduced in going to larger reactants, e.g. from +85 kJ mol−1 for B3LYP, Fig. 6, to +69 kJ mol−1 for B3LYP in Fig. 8.

On the other hand, the hydrogen abstraction in the intermediates, e.g. from [19] to [20] in Fig. 8, is connected with increasing barrier if systems become larger; as example we mention the barrier of +51 kJ mol−1 (B3LYP) relative to the educts from 5 to 6 as compared to +97 kJ mol−1 (B3LYP) from [19] to [20].

We attribute this increase to the fact that the reaction center is more crowded in the larger systems, which hinders the reaction.

We expect that both trends continue, even if reactions of large PAH are considered.


We have considered a variety of reactions to rationalize soot formation, e.g. under conditions typical for pyrolysis.

Our treatment is based on the assumption that radical formation, e.g. by hydrogen abstraction, is the initial step, on which the HACA mechanism is based too.

We did not consider the reaction PAH + C2H, since these require additional hydrogen migration and can only proceed by a perturbation of the aromatic system.

We proposed some new reaction pathways, Fig. 5, and steps for the HACA mechanism, Fig. 3 from [10] over [14] to [15] or from [11] over [12] and [13] to [15].

Our investigations confirm that the reactions, once a radical is available, proceed with relatively low barriers to make this good candidates for growth of PAHs and soot particles.

As a byproduct we confirm a previous result that PM3 is occasionally insufficient for a semi-quantitative description.

We have further considered condensation of two PAHs, one of them again assumed to be a radical.

For three typical reactions, Figs. 6–8, we have obtained reaction pathways, which proceed with barriers of typically less than 100 kJ mol−1.

This implies that these condensations could play an important role in soot formation.

The influence of entropies can easily be computed within the harmonic approximation; since we are dealing with radical reactions there could be a non-negligible contribution from the electronic term in the entropy.

We point out that the temperature dependency is in the ΔH term and the entropy is in the pre-factor of the Arrhenius equation, which in the modeling of kinetics is typically used as an adjustable parameter.

Lets come back to the temperature dependence mentioned in section 1, a temperature raise from 298 K to 1100 K has only a minor influence on ΔHr: about 6 kJ mol−1 for the first reaction of the acetylene addition and about 20 kJ mol−1 for the first aromatic condensation (benzene + naphthyl radical) (within the harmonic approximation).

We regard this again as not essential.

The calculations have been carried out at B3LYP level, which appears the standard choice for this type of investigations.

These results have been compared with those employing the BP86 functional, which always yields lower barriers but virtually identical reaction energies.

The deviations up to 35 kJ mol−1 in barriers and 28 kJ mol−1 in reaction energies appear not to be a serious problem under typical pyrolysis conditions (1300 K).

Since the non-hybrid functional requires typically 10% of CPU times necessary for B3LYP, this paves the way to treat systems with up to 100 atoms.