On the Importance of Clar Structures of Polybenzenoid Hydrocarbons as Revealed by the π-Contribution to the Electron Localization Function

Abstract

The degree of π-electron (de)localization and aromaticity of a series of polybenzenoid hydrocarbons (PBHs) has been analyzed through the π-contribution to the electron localization function (ELF_π), calculated at the B3LYP/6-311G(d,p) hybrid density functional theory level. The extent of π-electron delocalization in the various hexagons of a PBH was determined through analysis of the bifurcation values of the ELF_π basins (BV(ELF_π)), the spans in the bifurcation values in each hexagon (ΔBV(ELF_π)), and the ring-closure bifurcation values of the ELF_π (RCBV(ELF_π)). These computed results were compared to the qualitative description of local aromaticities of the different hexagons in terms of Clar structures with π-sextets. Benzene, [18]annulene, and thirty two PBHs were analyzed at their equilibrium geometries, and benzene and triphenylene were also analyzed at bond length distorted structures. In general, the description of PBHs in terms of Clar valence structures is supported by the ELF_π properties, although there are exceptions. For PBHs at their equilibrium geometries there is a clear sigmoidal relationship between the CC bond lengths and the amount of π-electron (de)localization at these bonds, however, this relationship is lost for bond distorted geometries. In the latter cases, we specifically examined benzene in D_3h symmetric “1,3,5-cyclohexatriene” structures and triphenylene in eight different structures. From the distorted benzenes and triphenylenes it becomes clear that there is a distinct tendency for the π-electron network to retain delocalization (aromaticity). The ELF_π analysis thus reveals an antidistortive rather than a distortive behavior of the π-electrons in these investigated compounds.

1. Introduction

In recent years, polybenzenoid hydrocarbons (PBHs) have received growing attention as a result of their increased number of applications in electronics and materials sciences. For example, pentacene and its derivatives are very promising materials for applications in organic thin-film transistors due to high field effect mobilities [1,2,3]. Graphene, isolated in 2004 as a single-layer sheet [1], has attracted booming interest from a range of different researchers due to its unique conductivity and hardness properties [4,5,6,7,8,9,10,11,12,13,14].

It is well known that the electronic structures of finite-size PBHs depend strongly on topology. For example, Clar observed in the 1950s that the relative stabilities and UV/Vis absorption maxima of a series of isomeric PBHs depend markedly on connectivity [15,16]. The heptabenzenoid isomers 1–6 (Figure 1) were examined and it was concluded that the isomer that could be written with a resonance structure having the largest number of π-sextets was the most stable and also had the most blue-shifted UV/Vis absorption (λ_β = 295 nm) [16]. Here, a π-sextet in a PBH is defined as a benzene ring with six delocalized π-electrons separated from adjacent benzene rings by formal CC single bonds, and the Clar structure is the valence structure which has the maximum number of disjoint π-sextets. In such valence structures, the CC bonds that are not part of π-sextets must be of unique CC double or single bond type [17,18]. Thus, among the heptabenzenoid isomers of Figure 1, tetrabenzanthracene (6) is the most stable isomer as it possesses the largest number (five) of aromatic π-sextets. In contrast, heptacene (1) is the least stable as it only can be drawn as containing one π-sextet. Moreover, this single π-sextet migrates between all seven hexagons, opposite to the five π-sextets of 6 which are drawn as confined to the five non-adjacent hexagons. The PBH 6 is an example of a fully benzenoid hydrocarbon since all C atoms are part of π-sextets.

Many different methods to assess the degree of aromaticity exist [19], and recently, Solà and co-workers concluded that aromaticity indices connected to electron delocalization properties are more revealing than those based on geometric, energetic, or magnetic properties [20,21]. Topological analyses of the electron localization function (ELF) [22], and in particular the π-contribution, ELF_π, can be used to directly link molecular electronic structure properties with (anti)aromaticity. This has been exploited by Santos et al., Malrieu et al., and us [23,24,25,26,27]. The ELF_π describes the kinetic energy destabilization at position r due to Pauli repulsion between same-spin π-electrons, and it takes values in the range 0 ≤ ELF_π(r) ≤ 1. The function is defined so that values close to one are found in regions where a π-electron pair is localized, whereas values close to zero indicate strong Pauli repulsion between same-spin electrons and are found in regions separating electron pairs.

A necessary condition for aromaticity is that the π-electrons are well delocalized over the allegedly aromatic circuit. By variation of the isovalue for ELF_π isosurfaces one can locate bifurcation points, i.e., the points at which two or several ELF_π basins merge (or split). The values of the ELF_π at these points are the bifurcation values (BV(ELF_π)s).

As measures of π-electron (de)localization, and thus, as aromaticity indicators, we now apply: (i) the BV(ELF_π)s at the individual CC bonds, (ii) the range in the bifurcation values of the ELF_π basins of a ring [ΔBV(ELF_π)], and (iii) the bifurcation value at which the ELF_π basins close to a single cyclic basin over the complete ring, i.e., the ring-closure bifurcation value [RCBV(ELF_π)]. As a consequence of their pronounced π-electron delocalization, aromatic rings should have small ΔBV(ELF_π) (the threshold used herein is ≤ 0.20), and their RCBV(ELF_π)s should be relatively high (herein ≥ 0.65). In contrast, antiaromatic annulenes at their equilibrium geometries have large ΔBV(ELF_π)s, revealing strong localization of π-electron pairs to distinct CC bonds, and their RCBV(ELF_π)s are very low (typically below 0.2). We earlier applied the ELF_π based indices to characterize a set of annulenes as aromatic or antiaromatic in their S₀ and T₁ states, respectively [26].

The character of the individual π-bonds are important, and we classify a CC π-bond as delocalized if its BV(ELF_π) is in the 0.65–0.96 range [28]. Only if all CC π-bonds of a hexagon are of delocalized nature do we categorize it as an aromatic π-sextet. If the BV(ELF_π) of an ELF_π bifurcation along a CC bond is larger than 0.96 or smaller than 0.65 we categorize this bond as a double or single bond, respectively, and we do not allow any such bonds in an aromatic π-sextet. Values somewhat smaller than 0.65 and larger than 0.96 are regularly found in π-conjugated polyenes having single- and double bond alternant structures [e.g., at B3LYP/6-311G(d,p) level planar s-cis-1,3-butadiene has BV(ELF_π)s of 0.998 and 0.399, respectively], and will also be found in non-aromatic, yet π-conjugated segments of the PBHs.

The local aromaticities of the individual rings in PBHs and the validity of Clar’s rule have previously been examined with a range of different quantum chemical techniques, including chemical graph theory [29,30,31,32], electron delocalization methods [33,34,35,36,37,38,39,40,41], ring current methods [42,43,44,45,46,47,48], nucleus independent chemical shifts (NICS) [34,49,50,51,52,53,54], and the harmonic oscillator model of aromaticity (HOMA) [33,34,55,56]. Based on these theoretical investigations it is generally accepted that Clar’s rule and Clar structures serve as useful rules-of-thumb for predicting stabilities of PBHs and variations in local aromaticities of different hexagons. Interestingly though, Balaban and Randic recently revealed that results from schemes taking into consideration so-called Fries resonance structures (Figure 2), i.e., those with the most number of double bonds shared between hexagons, even better agree with observed properties [31].

An interesting issue is whether the formation of Clar π-sextets is a result of the σ- or π-electron framework. Maksić and co-workers analyzed this using different partitioning schemes of the Hartree-Fock energies. They concluded that the driving force leading to Clar structures are favorable interactions in the σ-electron framework [57]. With regard to the π-electron framework in annulenes and polyenes it has been concluded based on thorough valence bond theoretical studies by Shaik, Hiberty and co-workers that it is distortive and prefers a CC bond alternant structure (in benzene D_3h symmetry) [58,59,60]. Yet, π-systems are simultaneously resonance stabilized, and for benzene this resonance stabilization is largest at the symmetric hexagon, although it has been concluded that for the π-network the distortive “Kekulé term always wins over the resonance term [60]. Hence, first when supported by the σ-framework does the hexagonal structure with the maximal resonance stabilization result.

In this regard, Angeli and Malrieu recently deduced that benzene displays both a distortive π-electron delocalization energy as well as an antidistortive cyclic delocalization energy [61,62], and at D_6h symmetry these energies are 16.35 and 1.45 eV, respectively. Upon distortion to a bond length alternant structure, there is a gain in the π-delocalization energy but it is outweighed by a reduction in the cyclic delocalization energy, allowing the preference of the σ-system for the perfect hexagon to prevail. Interestingly, the cyclic delocalization is still large in CC bond alternant benzene derivatives (“1,3,5-cyclohexatrienes”), as the one by Siegel and co-workers [63]. Indeed, such bond alternant benzenes retain a large portion of the cyclic delocalization energy of benzene, in line with earlier findings on conservation of aromaticity as observed in calculations of ring-currents, NICS, and energetic measures (e.g., hydrogenation energies) [64,65,66,67]. In particular, Corminboeuf, Schleyer, Mo and co-workers recently applied the block-localized wave function (BLW) scheme to compute NICS(0)_πzz, and found that a geometric distortion from D_6h benzene to a D_3h symmetric “1,3,5-cyclohexatriene” only reduced the NICS value from −36.3 to −28.2 [68]. However, effecting an “electronic” localization (blocking of π-conjugation) through the BLW approach led to a NICS(0)_πzz of zero. Using the ELF_π analysis it should now be interesting to examine the effects of geometry distortions of PBHs on the π-electron (de)localization and the degree of aromaticity of the various hexagons as determined by the ELF_π properties.

The first aim of the present study has been to probe what emphasis ELF_π gives to Clar representations of PBHs. Are there any variations in the π-(de)localization (aromaticity) between different π-sextets? What about relative stabilities of isomeric PBHs in terms of Clar structures? Are there compounds for which Clar’s rule fails? We have also used ELF_π to investigate bond length distorted PBHs. In such systems, will the π-electrons arrange similarly as in the equilibrium structure or will they adapt to the distorted geometry? Is the π-distortive effect manifested through the ELF_π properties, or will the π-delocalized and aromatic electronic structure be largely conserved, in line with a significant cyclic delocalization energy also at the bond distorted geometries?

2. Results and Discussion

As reference compounds for the subsequent analysis we first examined benzene, [18]annnulene, and a selection of oligophenyls where the latter compounds serve as references which can be represented by valence structures exclusively composed of π-sextets arranged in linear, branched, or cyclic chains. We then examined Clar’s rule by first investigating the heptabenzene isomers contained in the Clar and McCallum study of 1960 [16], and a series of gradually larger PBHs. Next, we explored the connection between the π-electron (de)localization as given by the BV(ELF_π)s and the CC bond lengths at the equilibrium geometries of the PBHs. We also examined the relationships between ΔBV(ELF_π) and different aromaticity indices. Finally, we investigated bond length distorted benzene and triphenylene, and noted how the π-electron distribution at the distorted geometries are described by Clar representations when based on an analysis of BV(ELF_π) and ΔBV(ELF_π). The CC bond lengths of all compounds are given in Supporting information Figure S1. All calculations were performed at the B3LYP/6-311G(d,p) hybrid density functional theory level.

2.1. Oligophenyls

The BV(ELF_π) and ΔBV(ELF_π) of D_6h symmetric benzene (7) serve as reference points for a perfect π-sextet, and at B3LYP/6-311G(d,p) level they are 0.911 and zero, respectively, similar to our previous values at OLYP/6-311G(d,p) level (0.906 and zero) [26]. Oligophenyls 8–12 constrained to planarity represent PBHs that are fully described by linear, branched, and cyclic chains of π-sextets having no leftover CC double bonds, and as such, they represent suitable compounds for probing purely the interaction between singly linked π-sextets. E.g., will meta- vs. para-connectivity in oligophenyls lead to differences in the (de)localization and interaction between π-sextets? With regard to sexiphenyl (12) and [18]annulene (13), these compounds are segments of several of the larger PBHs. Their features at D_6h symmetric structures are therefore interesting, although these structures are either a higher-order saddle point (12) or a conformer that is most likely not the lowest in energy (13) even though the relative energies of the conformers of 13 depend very much on the level of theory used [69].

As seen in Figure 3, biphenyl (8) has a nonzero ΔBV(ELF_π) and a rather small BV(ELF_π) for the single bond interconnecting the two phenyl groups, indicating a weak conjugation. The meta-connected terphenyl (10) has nearly identical ΔBV(ELF_π)s to that of 8 for both the terminal and central phenyl groups, whereas the para-connected terphenyl (9) has slightly larger ΔBV(ELF_π) of the central ring, revealing a tendency towards a quinoid valence structure. The star-shaped 1,3,5-triphenyl benzene (11) has a ΔBV(ELF_π) of the outer phenyl groups similar to that of biphenyl, whereas that of the central phenyl ring is zero. It is noteworthy that the BV(ELF_π) of the central ring is significantly lower than that of benzene, a feature also seen in larger PBHs (vide infra). With regard to the cyclic 12 the phenyl groups have a ΔBV(ELF_π) which resembles that of biphenyl. Thus, the conjugation between the benzene rings in oligophenyls leads to small, yet non-zero ΔBV(ELF_π)s. Finally, D_6h symmetric [18]annulene 13 has an even π-electron distribution indicated by a ΔBV(ELF_π) ≈ 0 revealing no π-electron pair localization for either type of CC bond. The BV(ELF_π)s are in general lower than in benzene, but are still well within our aromatic range.

2.2. Heptabenzenoids

In 1960, Clar and McCallum studied the stabilities and UV/Vis absorption spectra of a series of heptabenzenoids, these being heptacene (1), benzohexacene (2), two dibenzopentacenes (3 and 4), tribenzotetracene (5), and tetrabenzanthracene (6) in Figure 1 [16]. It was observed that the stability increased and the absorption bands shifted to shorter wavelength as the numbers of aromatic π-sextets were increased.

We now optimized 1–6 at (U)B3LYP/6-311G(d,p) level and the computed energies reproduce Clar’s experimental observations on the relative reactivities well as the heptabenzenoids become more stable the more π-sextets can be drawn (Figure 1). Indeed, there is a good linear correlation (R² = 0.949, see Figure S2 in Supporting Information) between the number of π-sextets and the relative energies of the various isomers corrected for the zero-point vibrational energy. An even better correlation (R² = 0.964) is found when the more recent M062X hybrid density functional method [70] is applied, although the magnitudes of the relative energies differ (see Table S1). This correlation clearly suggests that the conclusion of Clar and McCallum on the stabilities of 1–6 based on the different reactivities (kinetic stabilities) extends to the thermodynamic stabilities of these PBHs.

With B3LYP, the least stable isomer (1) is identified as a singlet biradical (<S²> = 0.76), consistent with the earlier study of Bendikov et al. [71], although it should be noted that a recent benchmark theoretical study, including a focal point analysis applied onto the results of a series of symmetry-restricted electron correlated calculations of gradually increased sophistication, revealed that 1 has a ¹A_g closed-shell singlet ground state [72]. Tetrabenzanthracene (6), having the maximum number of π-sextets (five) among the selected heptabenzenoids, is the isomer with the highest thermodynamic stability according to our calculations. The two isomers with three π-sextets (3 and 4) have similar relative energies. The calculated BV(ELF_π) and ΔBV(ELF_π) reveal that each of the hexagons of 2–5 (Figure 4) that can be drawn with localized π-sextets indeed fulfill our strict criteria for aromatic character, even though none of the rings have as low ΔBV(ELF_π) as the oligophenyls 8–12 (Figure 3). On the other hand, the outer CC bonds of the hexagons in 1–5 which share a migrating π-sextet have BV(ELF_π)s that represent delocalized π-bonds. However, these bonds do not form a closed loop and therefore do not constitute truly aromatic systems according to our criteria. For 1 it can be noted that there is a slight difference in the ELF_π properties between the open-shell and closed-shell solutions, but the trend is the same with the smallest ΔBV(ELF_π) in the central hexagon and the largest in the two outer ones.

The highest variation in BV(ELF_π) within a confined π-sextet (0.152) is found in 5, and although this value is significantly larger than that of benzene (0) and each of the oligophenyls (≤ 0.070), it is significantly smaller than that of naphthalene (0.237, vide infra). For 6, which can be described by one single valence structure with five confined π-sextets, the calculated BV(ELF_π)s and ΔBV(ELF_π)s are in agreement with its Clar structure. However, there are clear interactions across the CC “single” bonds that connect the π-sextets, and these interactions are stronger than in any of the oligophenyls. Interestingly, the two bifurcation points in 6 at which the full ELF_π basin in two steps separates into five 6π-electron basins occur at nearly the same ELF_π values (0.445 and 0.449; for the ELF_π isosurface at the last bifurcation point, see Figure 4). It is also noteworthy that the BV(ELF_π)s of the outer rings are closer to the value of benzene (0.911) than are the BV(ELF_π)s of the central ring, manifesting the similar observation in 1,3,5-triphenylbenzene 11. As the BV(ELF_π)s depend on the number of π-electrons in a ring [26], one could suspect a polarization of the π-electron distribution to the outer rings, but a natural population analysis of 6 shows that the π-electrons are essentially equally distributed between the five rings (5.9944 π-electrons in the inner ring and 6.0014 in each of the four outer rings). The overall smaller BV(ELF_π)s of the central ring should thus stem from its central position in this PBH, a feature which is even more pronounced in the larger PBHs.

2.3. Other Polybenzenoid Hydrocarbons

Various aspects of PBHs in terms of their Clar descriptions, e.g., extent of π-sextet confinement in fully benzenoid hydrocarbons and the effect on (de)localization of migrating π-sextets, were examined through the gradually larger compounds 14–34 (Figure 5 and Figure 6), gathered into different groups so as to facilitate comparisons.

Linear [n]acenes 14, 15, 17, 21, and 1: For naphthalene (14) it is known from experiments that the four C_αC_β bonds (1.376 Å [73]) are shorter than the other CC bonds, and our calculated C_αC_β bond length (1.374 Å) agrees with experiments. The BV(ELF_π)s reveal larger π-electron pair localization to these bonds than to any other CC bonds, and according to our criteria the four C_α-C_β bonds should be classified as CC double bonds, although their BV(ELF_π) is at the borderline of being classified as delocalized. Still, the ΔBV(ELF_π) of the hexagons of 14 (0.237) are moderately high and reveal a tendency to π-electron localization. Moreover, the ΔBV(ELF_π) of the outer 10π-electron loop is also substantial (0.221).

In anthracene (15), the two terminal hexagons display the largest ΔBV(ELF_π), and when compared to 14, the tendency for localization of the π-electrons into 1,3-butadiene segments in these rings is slightly more pronounced. Similar as in 1, the outer CC π-bonds of the central hexagon display delocalized BV(ELF_π), but these CC bonds are not linked into a cyclic loop, and consequently, we do not categorize it as a truly aromatic π-system. Moreover, the ΔBV(ELF_π) of the peripheral 14π-electron loop (0.304) has increased compared to 14, revealing some π-electron pair localization to certain bonds, these being the CC bonds in the 1,3-butadiene segments of the terminal hexagons. The features observed in 1 and 15 are naturally found in tetracene (17) and pentacene (21). Noteworthy, as the [n]acenes become longer the BV(ELF_π) of the CC bonds between hexagons become gradually lower, and it is thus questionable if any of the internal hexagons can be considered as aromatic. Instead, the ELF_π analysis reveals that the long linear [n]acenes are best described as two polyenes flanking each other, a description which is analogous to that put forward by Bendikov et al. [71]. In 1, the best π-delocalized circuit with the highest RCBV(ELF_π)s (0.634) and the near-lowest ΔBV(ELF_π) (0.345) is the 30π-electron periphery, and this is also the circuit with the strongest ring-current in linear [n]acenes according to Bultinck, Fowler and co-workers [47].

Earlier, it has been debated if the degree of aromaticity increases or decreases as one goes from the outer to the inner hexagon of a linear [n]acene [74]. Even though none of the hexagons can be considered as aromatic based on our criteria of the ELF_π properties, the ΔBV(ELF_π)s of both 1 and 21 decreases as one goes from the outer to the inner rings, supporting the view of a somewhat larger tendency towards aromaticity in the inner ring; a trend which was also observed in an earlier NICS study [53]. On the other hand, when based on the RCBV(ELF_π) the opposite trend is found for the hexagons of 1 as the BV(ELF_π)s of the inter-ring CC bonds become slightly lower as one moves to the central hexagon (Figure 4). For 1, it should also be noted that the ELF_π properties reveal no difference between closed-shell and open-shell solutions in the trends for the degree of aromaticity of the outer vs. inner hexagons, in contrast to what earlier has been concluded based on the NICS, HOMA and para-delocalization index (PDI) [75].

Anthracene (15) vs. phenanthrene (16): Phenanthrene (16) is 5.0 kcal/mol more stable than anthracene at B3LYP/6-311G(d,p) level, in line with earlier calculations (4.2 kcal/mol) at BLYP/TZ2P level [76]. However, Bader and co-workers used AIM analyses to question the explanation that the higher thermodynamic stability of 16 stems from better π-bonding in 16 than in 15 [77]. Instead, they argued that a non-bonded interaction between the two flanking H atoms in the bay region of 16, which is stabilizing rather than destabilizing, is the reason for the enhanced stability. However, Poater et al. recently disproved this argument using a quantitative energy decomposition analysis [76]. In line with the Clar structures our BV(ELF_π) and ΔBV(ELF_π) reveal a better π-electron delocalization in 16 than in 15 as two semi-confined π-sextets are formed in the A hexagons of 16 (Figure 6). Still, when comparing 16 with planar D_2h symmetric biphenyl (8), the ΔBV(ELF_π) of the hexagons in 8 is 0.052, lower than in 16. Interestingly, the lowest BV(ELF_π) in the A hexagons (the π-sextet rings) of 16 is found in the two inter-ring CC bonds between A and B hexagons, i.e., the two CC bonds that have double bonds in the Fries resonance structure. The difference in BV(ELF_π) to the two adjacent CC bonds in hexagons A is, however, minute.

Tetrabenzenoids 17–19: According to Clar’s rule triphenylene (19) as a fully benzenoid hydrocarbon should be the most stable of the three isomeric tetrabenzenoids 17–19, whereas tetracene (17) with one migrating π-sextet should be the least stable. Yet, the relative B3LYP/6-311G(d,p) energies are 9.5 kcal/mol for 17, -0.6 kcal/mol for 18, and 0.0 kcal/mol for 19. With regard to the π-electron (de)localization, the ΔBV(ELF_π)s reveal that 19 is well described by the Clar structure with three π-sextets, whereas chrysene (18) according to the ELF_π properties is most suitably described by a structure with two π-sextets confined to the two A hexagons. The π-electron structure described by the BV(ELF_π)s is in rather good accord with the Fries structure as the inter-ring CC bond between the B hexagons has a higher BV(ELF_π) than the adjacent bonds. A dominant description of 18 should be as 1,6-diphenyl-1,3,5-hexatriene.

Polybenzenoids 16, 18, 22–24: When going from phenanthrene (16) to chrysene (18), the π-sextets become slightly less confined to the two terminal rings, but when another hexagon is fused to chrysene, yielding picene (24), the π-electron delocalization in the terminal rings remains as in 18. However, the central ring of 24 has a large ΔBV(ELF_π) and a π-electron pair is localized to one of its CC bonds, inconsistent with its Clar structure according to which the central ring should also host a π-sextet. On the other hand, 1,2,5,6-dibenzanthracene (22) and 1,2,7,8-dibenzanthracene (23), which are almost isoenergetic with 24 (both 0.2 kcal/mol lower in energy), display three confined π-sextets. A lower aromaticity of the central ring of 24 than of 22 is in line with the observations of Solà and co-workers when based on PDI, NICS, and HOMA [34]. Clearly, there are differences in the π-electron (de)localization in the central rings of 22–24 which are not reflected by the Clar structures. Interestingly, there is no variation in the ΔBV(ELF_π)s of the central rings of 22 and 23 even though the three π-sextets have either para (22) or meta (23) connectivity, in disagreement with the central ring ΔBV(ELF_π)s of the two triphenyls 9 and 10 (Figure 3). A difference between 24 vs. 22 and 23 is that in 24 the C hexagon can share double bonds with both flanking B hexagons in the same resonance structure, and this Fries structure dominates according to ELF_π. Because of the connectivities of the C hexagons in 22 and 23, a similar double bond sharing as in 24 is not possible in these two PBHs.

Polybenzenoids 8, 16, 20, and 28: When going from phenanthrene (16) to pyrene (20), the BV(ELF_π) of the central CC bond reveals an increased π-electron delocalization between the two A rings, much stronger than in the idealized Clar structure with two π-sextets represented by biphenyl (8, Figure 3). Even though the BV(ELF_π) show a significant π-electron delocalization at this bond the two π-sextets are still rather confined to the A rings of 20. Indeed, recent studies of the Pauling-Wheland resonance energy (PWRE) per C atom showed that this energy gradually increases when going from benzene via biphenyl and phenantrene to pyrene from 10.2 to 12.3 kcal/mol, thus revealing the existence of an extra stabilizing non-Clar conjugation between the two idealized Clar π-sextets [78], in agreement with our BV(ELF_π) data. On the other hand, the ring-current map of 20 shows that the main current runs along the periphery [47], with no involvement of the central CC bond due to its symmetric position. Clearly, the most delocalized CC π-bonds are not necessarily parts of the strongest diatropic ring current paths, and in 20 the main ring current path even involves CC bonds with the largest variation in BV(ELF_π). Pyrene is thus an interesting compound for which aromaticity indices based on the electron density as compared to those based on magnetic properties give significantly different answers as to which cycles are the aromatic ones.

Next, when the two localized CC double bonds of the B rings become part of hexagons in the fully benzenoid hydrocarbon 28, the confinement of all four π-sextets increases considerably. Interestingly, the BV(ELF_π) of the internal CC bond is also reduced when compared to that of 20. Clearly, the ELF_π analysis supports that fully benzenoid hydrocarbons are well described by Clar representations. The ΔBV(ELF_π), however, reveals a difference between the A and C hexagons with the ΔBV(ELF_π) of the A hexagon being slightly smaller.

Polybenzenoids 16, 19, 27, and 28: A significant increase in π-sextet localization is also found when going from 16 to triphenylene (19). However, when an ethylene moiety is added to 19, the π-electron delocalization of the A hexagons in 27 is reduced so that the A, B, and C rings form a segment with similar ELF_π properties as 20. The strongly confined π-sextet character is regained in 28.

Polybenzenoids 14, 25, 26, 29, and 30: Perylene (25) shows a clear resemblance with naphthalene (14) as the coupling between the two segments is weak (Figure 6), and this observation is in line with the earlier calculated ring currents [47]. The next species, dibenzoperylene (26), shows a completely different delocalization pattern which can be viewed as a combination of the BV(ELF_π) patterns of pyrene 20 and the para-connected terphenyl 9, although the central C hexagon of 26 reveals no tendency towards a quinoid resonance structure. Interestingly, the pattern for π-electron delocalization according to ELF_π differs significantly from the dominant ring-current path which is found in a circuit around the periphery [47,48]. It is particularly notable that the major ring-current path involves four CC equivalent bonds with BV(ELF_π)s of merely 0.525, revealing once more that the predominant ring-current path does not necessarily coincide with the best π-electron delocalized path.

When an ethylene moiety is added to 25, yielding 29, π-sextets become well-confined to the two C hexagons, however, with regard to the third π-sextet which can be drawn in 29 it is only semi-confined to hexagon B even though all CC bonds of this hexagon have BV(ELF_π)s within the π-delocalized range. To draw 29 as a Clar structure with three π-sextets in accordance is thus a borderline case according to the ELF_π results.

When one more ethylene moiety is added to 29, coronene (30) is formed, a PBH which could be described by three delocalized π-sextets. Initially, 30 was visualized as superaromatic with an [18]annulene circumscribing a benzene ring [18]. However, later calculations of the current-density map showed opposing ring currents; a diatropic ring current in the 18-carbon rim and a paratropic ring current in the 6-carbon hub [79,80]. Now, with regard to the π-electron (de)localization the BV(ELF_π)s in the 18-carbon rim (0.967 and 0.713) are far from those of D_6h symmetric [18]annulene (0.867 and 0.878, Figure 3), and they reveal a clear tendency for π-electron pair localization to some of the peripheral CC bonds. This also agrees with the previous MP2 results of Ciesielski et al. [80].

Polybenzenoids 31–34: The hexabenzocoronene isomer 31 is a fully benzenoid hydrocarbon with all C atoms included in π-sextets, and the ELF_π results support this view. Also the Clar structures of circumcoronene (32) and kekulene (33) are supported by BV(ELF_π)s and ΔBV(ELF_π)s as the A and D (32) and A (33) hexagons host π-sextets. Kekulene, synthesized by Staab and Diederich in 1978 [81], was earlier suggested to be a benzenoid rather than a superaromatic annulenoid hydrocarbon according to graph theory, NICS, and current-density map calculations [82,83,84]. From the ELF_π results one can note that both the inner [18]annulene and the outer [30]annulene paths have large variations in BV(ELF_π) suggesting poor overall π-delocalization along these circuits. With regard to 33, the idealized BV(ELF_π) and ΔBV(ELF_π) now suggest that the π-sextets in this PBH are less isolated, and thus somewhat less aromatic, than those in 12, the idealized Clar structure for six π-sextets in a cycle.

If one regards 7, 11, 31 and 34, i.e., the fully benzenoid hydrocarbons, two interesting items can be observed. First, if one compares the BV(ELF_π)s of the most central and second most central CC bonds of 31 and 34 one finds that the difference becomes smaller (0.179 in 31 and 0.106 in 34, respectively). Naturally, the difference between the various hexagons in a large PBH will eventually disappear, however, this difference seems to be diminished at the cores of rather small PBHs. The difference between the hexagons at the core of PBHs with 40–50 hexagons is likely very small. Secondly, the BV(ELF_π) of the central hexagon gradually decreases from 0.911 in 7 to 0.741 in 34. For the uniform electron gas, the ELF_π should be everywhere 0.5 [22], and even though perfect π-conjugation will not lead to a uniform electron gas, the larger and the more perfectly conjugated a π-system becomes, the lower the BV(ELF_π)s seem to become.

2.4. Comparison with Previous Aromaticity Indices

Many of the PBHs investigated herein were also examined earlier by various groups using a series of different aromaticity indices, and, although ΔBV(ELF_π) is not suitable as a general aromaticity index because it depends on the ring size [26], the data presented herein allows for an assessment of ΔBV(ELF_π) against these earlier indices. For the analysis we compared against the cyclic resonance energy (CRE) of Aihara [39], the Atoms-in-Molecules based properties of the electron density at the ring critical points (here λ_π,3) of Howard and Krygowski [31], as well as the bond order index of aromaticity (BOIA), six-center bond index (SCI), HOMA and NICS(0) data reported by Bultinck and co-workers [85] (see Table S2, Supporting Information). Note that our comparisons are against NICS(0) data, i.e., the very first 1996 version of NICS according to which the isotropic NICS value is computed at the ring center of a cyclic conjugated compound [49]. Since then, significant NICS refinements have been introduced [50,52,54]. The most sophisticated NICS(0)_πzz version, based on only the contributions of the zz (perpendicular) tensor components of the π-MO’s, gives the best results, although the quality of the more readily available NICS(1)_zz data (1 Å above the ring center) also is very high.

Among the indices to which we compare, the BOIA, SCI and λ_π,3 are directly linked to properties of the electronic structure, NICS is a magnetic index, HOMA is a geometric, and CRE is a magnetic-energetic based index. Significant correlations between ΔBV(ELF_π), BOIA, SCI, and λ_π,3, as electronic structure based aromaticity indices, can be expected. The comparisons are based on eighteen PBHs for BOIA, SCI, NICS(0) and HOMA, and ten PBHs for λ_π,3 and CRE. A comprehensive comparison between several of the previously reported indices has been presented by Bultinck [86]. Moreover, the connection between classical (structural and energetic) and magnetic aromaticity indices has also been analyzed earlier, and it has been concluded that there are no discrepancies between the various classes of measures provided that one compares values from indices which report on the truly local aromatic character of the hexagons in a PBH [87,88].

The relationships between our ΔBV(ELF_π) data and the data from the six earlier aromaticity indices are displayed in Figure 7. Interestingly, in each of the graphs one can note a few very distinct outliers (hexagons type II), marked in black. These outlying data points correspond to hexagons that do not host Clar π-sextets and they are found in 6, 19, 25, 27, 28, and 30, i.e., small and symmetric PBHs. The hexagons type II include the non-Clar hexagon in 6, hexagon A in 19, hexagon B in 25, hexagon C in 27, hexagon B in 28, and hexagon B in 30 (Figure 5). Thus, the hexagons of type II are in all except one case (hexagon B in 30) situated between strongly aromatic segments that host Clar sextets, and therefore, one could expect them to display low local aromaticity. This is also the situation with ΔBV(ELF_π), and indeed, BOIA, SCI, λ_π,3, and CRE describe them as even less aromatic than ΔBV(ELF_π). On the other hand, HOMA describes them as significantly more aromatic, nearly as aromatic as benzene, so that HOMA comes at conflict with both the Clar description and all aromaticity indices based on electronic structure. The “aromatic” HOMA values for hexagons of type II likely stem from geometric constraints rather than from true aromaticity.

However, for the clear majority of the data points (hexagons type I) there are good correlations between ΔBV(ELF_π) and λ_π,3, BOIA and HOMA (R² ≥ 0.90), as well as moderate correlations with SCI and CRE (0.70 ≤ R² ≤ 0.90). These correlations are of similar quality to those found between most of the indices recently examined by Bultinck [86]. On the other hand, there is no correlation between ΔBV(ELF_π) and NICS(0), however, a few points of concerns should be raised. Clearly, the comparison of the ELF_π data should be made against the more recent and refined NICS versions, i.e., NICS(0)_πzz or NICS(1)_zz [54], rather than against data from the old NICS(0). Recent studies have concluded that NICS does not reveal the local aromaticity of individual rings of a PBH but rather that of several circuits combined [46,47,89,90], although there should be no inconsistencies between NICS and local aromaticity measures provided that interfering ring-current contributions of contaminating circuits are taken care of [88]. However, to further expand the understanding of NICS an assessment of the performances of the various recent NICS refinements for the determination of the local aromaticity in PBHs should also be valuable. In the context of magnetic vs. structural/energetic aromaticity measures, it should be noted that in addition to the poor correlation between ΔBV(ELF_π) and NICS(0), the correlation between the ΔBV(ELF_π) and CRE data is slightly weaker than those between ΔBV(ELF_π) and the remaining indices.

2.5. The ELF_π Bifurcation Values and CC Bond Lengths

For the PBHs 1–7 and 14–30 at their equilibrium geometries, we plotted the BV(ELF_π)s against the CC bond lengths and fitted a sigma function to the data set (Figure 8). As the BV(ELF_π)s range from 0 to 1, and as these two end values are approached in regions between two atoms which either separates two electron pairs (BV(ELF_π) → 0) or host an electron pair (BV(ELF_π) → 1), we selected the sigma function for the fit of BV(ELF_π)s against CC bond lengths, although there is no strict justification for this selection. PBHs 8–12 are model compounds constrained to planarity, [18]annulene (13) is not a PBH, and PBHs 31–34 were found to have minute but existing deviations from planarity (see Figure S3), and therefore these compounds were excluded from the function fit. Although there is no strict justification for a sigmoidal dependence of the BV(ELF_π)s on the CC bond lengths, there is a reasonable agreement between the fitted function and the calculated values in the bond length interval 1.350 – 1.476 Å covered by the 23 PBHs [91]. This indicates that there is a direct relationship between the extent of π-electron (de)localization at each of the CC bonds, as given by the BV(ELF_π)s of the corresponding bifurcations, and the optimal CC bond lengths.

In order to test if the sigmoidal dependence of the BV(ELF_π)s on CC bond lengths is valid also for [4n]annulenes, and to extend the bond length range, we insert the data points for D_2h and D_4h symmetric cyclobutadiene (CBD) and cyclooctatetraene (COT), respectively, into Figure 8. Three of the four data points fall perfectly onto the sigmoidal curve, but there is a slight exaggeration of the BV(ELF_π) for the longest CC bond of CBD.

Shaik, Hiberty and co-workers, have revealed that forces within the σ-framework of benzene lead to the observed D_6h symmetric structure, whereas the six π-electrons exert a π-distortive effect that strives towards a bond alternant D_3h symmetry [58,59,60]. Furthermore, Maksić and co-workers deduced that σ-type interactions are also the major contributors to the Clar structures of PBHs, and they showed that these structures are the only minima on the potential energy surfaces [57]. Consequently, if the π-distortive effect dictates the behavior of the π-electron distribution so that it prefers a bond alternant structure, then the π-electron distribution should readily adjust if the geometry of a PBH is made more bond alternant than at its equilibrium geometry. If this is the case, the BV(ELF_π)s of the bifurcations at the CC bonds of distorted PBHs should follow the sigmoidal curve.

We now examined benzenes with D_3h symmetry (the b_2u distortion mode) for which the CC bonds were set to alternate at 1.40 ± 0.01 Å until 1.40 ± 0.08 Å, respectively. To probe the effect of geometry distortions of a PBH, we also examined the ELF_π properties of the fully benzenoid hydrocarbon 19 set at the bond length distorted geometries 19-II–19-IX. Triphenylene was chosen for this analysis as it is a small fully benzenoid hydrocarbon which is well described by its Clar structure according to the ELF_π analysis. Among the distorted triphenylenes, 19-II has every CC bond fixed at 1.400 Å, 19-III and 19-IV have the three internal inter-ring CC bonds extended to 1.500 and 1.600 Å, respectively, and 19-V and 19-VI have all CC bonds of the central hexagon set at 1.400 Å, whereas in the three outer hexagons the bond lengths are 1.500 and 1.300 Å (19-V) and 1.600 and 1.200 Å (19-VI), respectively. We also probed three triphenylenes 19-VII–19-IX were all CC bond lengths alternate at 1.350 and 1.450 Å. Two of these (19-VII and 19-VIII) allow for π-sextet formation in the three outer rings, whereas in the last (19-IX) the bond lengths are set so that the bond lengths match the anti-Fries structure (Figure 2), and as a consequence, the π-sextet should form in the central hexagon.

Interestingly, for the distorted benzenes (purple marks, Figure 8), the BV(ELF_π)s at CC bonds longer than 1.400 Å correspond to more delocalized π-electron distributions than anticipated based on the sigmoidal curve. On the other hand, the BV(ELF_π)s of bifurcations at bonds shorter than 1.400 Å correspond to less distinctly formed π-electron pairs than is ideal for this CC distance according to the sigmoidal curve. Thus, the ELF_π properties reveal that the π-electron network adheres to a delocalized (aromatic) arrangement even at geometries that correspond to a bond length alternant “1,3,5-cyclohexatriene”. Noteworthy, for the D_3h symmetric structure with CC bonds of 1.32 and 1.48 Å, the ΔBV(ELF_π) equals 0.293, much smaller than the value expected from the sigmoidal fit (0.722). We even examined benzene at an extreme “1,3,5-cyclohexatriene” geometry with CC “double” and “single” bonds of 1.200 and 1.600 Å, respectively, giving a ΔBV(ELF_π) of merely 0.730, in contrast to the value of ~1 estimated from the fitted function. According to our aromaticity criterium ΔBV(ELF_π) ≤ 0.2, benzene becomes “nonaromatic” at a bond alternant geometry with bond lengths of approximately 1.455 and 1.345 Å. These findings clearly fall in line with the observed strong calculated ring currents for constricted benzenes [64,65,67], such as the one of Siegel and co-workers [63], and with the NICS and energetic data of Stanger as well as Corminboeuf and colleagues [66,68]. Consequently, even when given the opportunity through the help of geometric distortions, an aromatic π-electron system is not prone to distort, but in contrast adheres to a delocalized distribution which is not very different from that of the maximally aromatic equilibrium geometry.

A similar tendency to adhere to a π-electron delocalized arrangement in the lowest ππ* excited triplet state of COT was earlier found by us through the ELF_π properties [26], in line with the triplet state aromatic character of 4n π-electron annulenes [90,92,93,94,95,96,97]. In its lowest ππ* triplet state (T₁), COT at OLYP/6-311G(d,p) level has D_8h symmetry, a ΔBV(ELF_π) of zero, and a RCBV(ELF_π) of 0.845 [26]. When the CC bonds of COT in the T₁ state were constrained at the optimal bond lengths of the S₀ state structure with D_4h symmetry (1.353 and 1.475 Å, respectively), a ΔBV(ELF_π) of only 0.116 was calculated, whereas in the anti/nonaromatic S₀ state the ΔBV(ELF_π) was as large as 0.637. This is clearly analogous to what we found above for benzene in the electronic ground state.

To gain insight into if and how π-electron interactions influence the Clar structure of a PBH, we examined the BV(ELF_π)) and ΔBV(ELF_π) for both optimized and bond-length distorted triphenylene (Figure 9 and Figure 10). The equilibrium geometry (19-I) is D_3h symmetric and all distorted structures (19-II–19-VI) keep this symmetry. When all CC bonds are set at 1.400 Å (19-II), the π-electron system is still very well described by the original Clar structure with three π-sextets. As a consequence of the Clar type π-electron structure the BV(ELF_π)s of the CC bonds of 1.400 Å vary significantly, and as seen in Figure 9, they are often far from the ideal value as given by the sigmoidal curve. In 19-III, all CC bonds of the three outer hexagons are still within the delocalized range, even though ΔBV(ELF_π) is larger than the 0.2 used herein as a threshold for aromaticity. When the internal CC bonds are expanded to 1.600 Å (19-IV), the triphenylene tends towards the character of an [18]annulene rather than of a triphenylene with three π-sextets, even though the BV(ELF_π) value at the internal CC bonds of 1.600 Å is as high as 0.512, far from the ideal value close to zero. It is thus possible to disrupt the Clar structure arrangement of the π-electrons in 19 to some extent, but only through extensive distortions. The relative energies of 19-III and 19-IV are 12.1 and 31.4 kcal/mol, respectively.

In the second type of distortions, represented by 19-V and 19-VI, the central hexagon is constrained to aromatic CC bond lengths of 1.400 Å, and it is seen that in both structures this hexagon can be classified as aromatic. Despite the central aromatic hexagon the CC bonds of the outer hexagons still adhere to a more π-delocalized situation than ideal for the particular CC distances. Noteworthy, 19-V and 19-VI are energetically very unfavorable as they are 38.5 and 182.8 kcal/mol above 19-I. Among the last three distorted triphenylenes (19-VII–19-IX), the Clar structural description is maintained in structures 19-VII and 19-VIII whereas 19-IX with short outer CC bonds of the central hexagon, a π-sextet has formed in the central ring, even though there is still a significant delocalization in the three outer rings. Also for 19-VII–19-IX, the BV(ELF_π)s clearly indicate that the aromaticity is retained as their values are far from the sigmoidal curve.

Consequently, the π-electrons in PBHs seem to “lead a life” which, at least to some extent, does not conflict with the desire of the σ-framework for Clar structures observed by Maksić and colleagues [57]. The π-electrons are willing to abandon their clustering into π-sextets, but only after brutal force from the outside has been applied as these distorted geometries are of high relative energies. Similar as there is an antidistortive cyclic delocalization energy in benzene which is maximal at the hexagonal structure [61], there is likely an extended cyclic delocalization energy in PBHs which favors formation of π-sextets and Clar structures. Taken together, the ELF_π data presented here indicate that there is a π-antidistortive effect which is operative in the PBHs.

3. Computational Methods

Geometry optimizations were performed using the B3LYP hybrid density functional theory (DFT) method [98,99], together with the 6-311G(d,p) valence triple-zeta basis set of Pople and co-workers [100]. Frequency calculations were performed to check the nature of the optimized structures (minima or transition states). The relative energies of the heptabenzenoids 1–6 were also examined at the recently developed M062X hybrid density functional theory level of Zhao and Truhlar [70] in order to place the corresponding B3LYP energies on a more solid foundation. All electronic structure calculations were carried out using the Gaussian 03 package [101].

The electron localization function (ELF) was introduced by Becke and Edgecombe as a tool for topological analysis of atomic and molecular electronic structure [22], and it is expressed as:

$$ELF(r)={\left[1+{({\chi }_{\sigma }(r))}^2\right]}^{-1}={\left[1+{\left[\frac{T_{\sigma }(r)}{T_{\sigma }^0(r)}\right]}^2\right]}^{-1}$$

where χ_σ is a dimensionless localization index by which the excess local kinetic energy, T_σ(r), due to Pauli repulsion between σ-electrons (σ = α or β) at position r is calibrated against the Thomas-Fermi kinetic energy of a uniform electron gas. With this definition, the ELF takes values in the range 0 ≤ ELF(r) ≤ 1, where values close (equal) to one are found in regions of space with strong (perfect) electron localization, and values close to zero are found in regions with large excess kinetic energy such as in boundary regions between two electron pairs where same-spin electrons come close together. The σ- and π-components of the kinetic energy can be separated, and this also applies to the ELF, even though the sum of the two components (ELF_σ and ELF_π) does not correspond to the total ELF. The ELFs, including only the orbitals involved in the π-system, giving ELF_π, were calculated with the Dgrid program [102], and the ELF isosurfaces were visualized with the VMD package [103].

The compounds 1–7, 13–30 were planar at their equilibrium structures, whereas 8–12 are model compounds constrained to planarity, and 31–34 had minute but existing deviations from planarity at their equilibrium geometries. Still, the analysis of 31–34 was performed using the planar structure.

The CC bond lengths and Cartesian coordinates of the optimized structures, absolute energies, symmetries, and number of imaginary frequencies at B3LYP/6-311G(d,p) level for 1–34 are provided in the supplementary material.

4. Conclusions

The degree of aromaticity of benzene, [18]annulene, and a series of polybenzenoid hydrocarbons (PBHs) (1–34) has been analyzed through the π-contribution to the electron localization function (ELF_π) at the B3LYP hybrid density functional theory level. The ELF_π results show that for most PBHs which can be described with a Clar structure having all π-sextets confined to single hexagons, their Clar structures are in general reproduced well, and this applies in particular to the fully benzenoid hydrocarbons. For the PBHs with migrating π-sextets, e.g., the oligoacenes, the ELF_π analyses reveal that the outer CC bonds in the hexagons with the migrating π-sextets have delocalized π-bonds, but that none of these rings separately are truly aromatic as the inter-ring CC bonds have low BV(ELF_π)s.

It should however be noted that the cyclic pathways with the smallest variation in BV(ELF_π) are not always the paths that host the largest diatropic ring currents. E.g., pyrene (20) displays the main current in the peripheral CC bonds, although these bonds display significant variation in the π-(de)localization according to ELF_π. Clearly, ring-current based indices will not always provide a similar gradation of local aromaticity of the rings in a PBH as indices based on properties of the electron density. In this context, the performance of the recent NICS refinements such as NICS(1)_zz and NICS(0)_πzz for the determination of local aromaticities in PBHs should be assessed. Another finding of the present study concerns the relative energies of isomeric PBHs. For heptabenzenoids 1–6 the calculated relative energies clearly correlate with the number of π-sextets, however, such a relationship is less apparent for the smaller 17–19 and 22–24.

For the equilibrium geometries there is a sigmoidal relationship between the CC bond BV(ELF_π)s and the corresponding CC bond lengths. Based on this sigmoidal relationship we examined bond length distorted structures. From an analysis of bond distorted benzene and triphenylene we find that the π-electron system of aromatic sextets retain a more delocalized distribution than what is normal for the particular bond lengths, i.e., the π-electron system of an aromatic π-sextet displays antidistortive behavior. Indeed, if the π-system would be distortive and prefer a symmetry which allows for bond length alternation, then the BV(ELF_π)s of the distorted geometries should follow the sigmoidal curve. However, the large divergences for the distorted triphenylene structures suggest something different; the π-electron system also prefers to arrange in Clar structures.