Aromaticity Determines the Relative Stability of Kinked vs. Straight Topologies in Polycyclic Aromatic Hydrocarbons

It is well-known that kinked phenacenes are more stable than their isomeric linear acenes, the archetypal example being phenanthrene that is more stable than anthracene by about 4–8 kcal/mol. In previous studies, the origin of the higher stability of kinked polycyclic aromatic hydrocarbons (PAHs) was found to be better π-bonding interactions, i.e., larger aromaticity, in kinked as compared to linear PAHs. Some years ago, however, Dominikowska and Palusiak (2011) found that dicationic linear anthracene is more stable than the dicationic kinked phenanthrene. Therefore, these authors showed that, in some cases, the linear topology in PAHs can be preferred over the kinked one. Our results using energy decomposition analyses in combination with the turn-upside-down approach show that the origin of the higher stability of dicationic anthracene is the same as in the neutral species, i.e., better π-bonding interactions. A similar result is found for the kinked and straight pyrano-chromenes. We conclude that the aromaticity is the driving force that determines the relative stability of kinked vs. straight topologies in PAHs.


INTRODUCTION
The change in topology when going from linear anthracene to its isomeric kinked phenanthrene has important consequences. From an energetic point of view, phenanthrene is more stable than anthracene by about 4-8 kcal/mol (Balaban, 1980;Paniagua, 1991, 1993;Behrens et al., 1994;Kato et al., 2002;Matta et al., 2003;Randić, 2003;Poater et al., 2007b). As compared to anthracene, phenanthrene has a higher ionization potential (by about 0.4 eV) (Boschi et al., 1974;Dabestani and Ivanov, 1999;Kato et al., 2002) and a reduced electron affinity (by ca. 0.5 eV) (Tschurl et al., 2006), as well as, a larger HOMO-LUMO gap as indicated by theoretical calculations (Kato et al., 2002;Poater et al., 2007b). Moreover, there is a blue shift in the S 1 → S 0 transition when going from anthracene to phenanthrene (Dabestani and Ivanov, 1999). All these observations indicate a higher kinetic and thermodynamic stability of phenanthrene with respect to anthracene. Still, both isomers undergo addition reactions usually at 9,10 position (Clar, 1964;Biermann and Schmidt, 1980;Wiberg, 1999). Theoretical calculations show that the ring currents are mainly localized in the central six-membered ring (6-MR) in anthracene, whereas for phenanthrene they are more intense in the outer 6-MRs (Anusooya et al., 1998;Ligabue et al., 1999;Steiner et al., 2002). These different ring currents lead to different NMR response and magnetizability values (Ligabue et al., 1999).
The origin of the larger relative stability of phenanthrene as compared to anthracene is controversial. Several authors (Matta et al., 2003;Bader, 2006;Wolstenholme and Cameron, 2006;Hernández-Trujillo and Matta, 2007;Eskandari and Van Alsenoy, 2014;Monteiro and Firme, 2014) based on the results from the quantum theory of atoms in molecules (QTAIM) consider that the H···H attraction between the hydrogen atoms of the bay region of phenanthrene (H4 and H5 attached to C4 and C5 in Scheme 1) justifies its higher stability. These two hydrogen atoms are about 5 kcal/mol stabilized with respect to hydrogen atoms in anthracene. This result was interpreted considering the presence of a H···H bonding interaction that contributes to the stabilization of phenanthrene by ca. 10 kcal/mol. Presence of this H···H stabilizing interaction would explain the larger stability of kinked phenacenes as compared to their isomeric linear acenes. It is worth mentioning that this interpretation differs from that of other authors that contemplate the existence of London dispersion stabilizations in hydrocarbons due to CH···HC weak interactions (Echeverría et al., 2011(Echeverría et al., , 2017Danovich et al., 2013).
Another more classical interpretation of the higher stability of phenanthrene derives from the analysis of the π-bonding interactions. Already in 1933, Pauling and Sherman found that the resonance of the π-system is more effective in phenanthrene than in anthracene (Pauling and Sherman, 1933). Later on, Dewar and De Llano (1969) as well as, Hess and Schaad (1971) confirmed these results by finding large resonance energies in phenanthrene as compared to anthracene (Cyranski, 2005). Fukui also explained the larger stability of phenanthrene by constructing the two benzenoids from naphthalene and butadiene fragments and showing that π-interactions are more stabilizing in phenanthrene (Fukui, 1982a,b). Clar's π-sextet model also indicates a larger stability of the π-system of phenanthrene by assigning two πsextets to phenanthrene and only one migrating π-sextet to anthracene (see Scheme 1; Clar, 1972;Randić, 2003;Portella et al., 2005;Solà, 2013). More recently, some of us (Poater et al., 2007b), using energy decomposition analyses proved that π-electron interactions between two 2-methtriyl-phenyl fragments (Ant or Phe in Scheme 2) that form phenanthrene are more stabilizing than those generated in the formation of anthracene using the same fragments. Moreover, by comparing the energies of the triplet states 4,5-didehydrophenanthrene and 4,5-didehydroanthracene (Scheme 1), we also found that the 4,5-didehydrophenanthrene species is 5.2 kcal/mol more stable than 4,5-didehydroanthracene (Poater et al., 2007b), whereas at the same level of theory phenanthrene is only 4.2 kcal/mol more stable than anthracene. Moreover, comparison between different didehydrophenanthrenes in their triplet states indicate that 4,5-didehydrophenanthrene is particularly stable (Poater et al., 2007a). In our opinion, these results are a clear indication of the repulsive character of the H···H interaction between the bay H atoms in phenanthrene. Thus, our work confirmed the classical picture of better π-bonding causing phenanthrene's enhanced stability, despite the unfavorable H···H repulsion (Poater et al., 2006a,b) between the hydrogen atoms of the bay region. This interpretation was reinforced by means of natural bond orbital (NBO) (Alkorta et al., 2008;Weinhold et al., 2014) and quantum kinetic energy density analyses (Jacobsen, 2010). Remarkably, Grimme et al. (2009) provided conclusive experimental evidence supporting the traditional view of steric repulsion between the H atoms of the bay region of phenanthrene. Finally, in 2007, Martín Pendás et al. (2007 in an attempt to reconcile the orthodox QTAIM interpretation with the classical view, showed that bond paths in QTAIM are not necessarily a sign of attractive or repulsive interactions, but indications of the presence of preferred quantum-mechanical exchange channels. In a more recent work, Dominikowska and Palusiak (2011), found that the anthracene dication is more stable than phenanthrene dication by about 16 kcal/mol. Therefore, by just removing two electrons from anthracene and phenanthrene, the relative stability of these two isomers is reversed. This result cannot be predicted in the framework of the conventional QTAIM H···H attractive interaction because this interaction should not be affected significantly by removing two π-electrons. On the other hand, it is easily understood by Clar's πsextet model and the classical picture of H···H repulsion. Both anthracene and phenanthrene dications have two π-sextets (see Scheme 1). From the Clar model, the stabilization of the πsystem is expected to be similar for the two C 14 H 2+ 10 isomers and, in this situation, the H···H steric repulsion should be the factor that makes phenanthrene dication less stable than anthracene dication. However, it is very important to take into account other factors not included in the Clar model to explain the phenanthrene dication reduction of stability. In particular, the loss of aromaticity due to partial delocalization of the positive charge to the external rings of phenanthrene plays also a very important role in explaining this reduction of stability (vide infra).
The example of the change in relative stabilities when going from neutral to dicationic phenanthrene and anthracene is an indication that the stability of the π-system, i.e., the aromaticity, plays a key role in determining the relative stability of kinked vs. straight topologies in polycyclic benzenoids. The main goal of this paper is to test this hypothesis by studying first the isomerization energies of the couples phenanthrene/anthracene and dicationic phenanthrene/anthracene, shown in Scheme 2. And next two couples of pyrano-chromenes: tetramethylated (linear 2,2,8,8-tetrametyl-2H,8H-pyrano[3,2-g]chromenehereafter anthraceneO2-vs. kinked 3, 3,8,8-tetramethyl-3,8dihydropyrano[3,2-f]chromene-hereafter phenanthreneO2-) and dicationic pyrilium derived (linear pyrano[3,2-g]chromene-1,9-diium-hereafter anthracenePyr-vs. kinked pyrano[3,2f]chromene-4,7-diium-hereafter phenanthrenePyr-), shown in Scheme 3. In the first pair of isomers, phenanthrene is more aromatic and more stable than anthracene, whereas in the second and third couples the aromaticity of the two isomers is expected to be similar because the number of π-sextets is the same for the two isomers and, therefore, it is likely that the straight isomer will be the most stable. For the last couple, we expect a behavior closer to that of anthracene/phenanthrene. Our investigation will be SCHEME 1 | Clar's representation of neutral, 4,5-didehydro, and dicationic anthracene and phenanthrene, and straight and kinked tetramethylated-pyrano-chromenes. The neutral representation of anthracene and phenanthrene with the labeling of the different carbon atoms.
performed by means of an energy decomposition analysis of the Kohn-Sham density functional theory (DFT) wavefunction using the turn-upside-down approach (El-Hamdi et al., 2011, 2013 to analyze the studied isomerization energies. We anticipate here that our calculations show that the aromaticity of straight and kinked isomers is the main factor that determines the relative stability of kinked vs. straight topologies in PAHs.

MATERIALS AND METHODS
All geometry optimizations were performed with the Gaussian 09 package (Frisch et al., 2009) by using the B3LYP hybrid density functional (Lee et al., 1988;Becke, 1993) and the 6-311++G(3df,3pd) basis set (Frisch et al., 1984) without symmetry constraints (Table S1). Analytical Hessians were computed to confirm that the optimized structures are indeed minima (zero imaginary frequencies).
Aromaticity was evaluated at the same level of theory by means of electronic-based aromaticity criteria Feixas et al., 2015) such as the multicenter electron sharing indices (MCI) (Giambiagi et al., 2000;Bultinck et al., 2005), the para delocalization index (PDIs) (Poater et al., 2003), and the aromatic fluctuation index (FLU) (Matito et al., 2005), as well as, with the geometry-based harmonic oscillator model of aromaticity (HOMA) (Kruszewski and Krygowski, 1972;Krygowski, 1993) and the magnetic-based nucleus-independent chemical shift (NICS) obtained using the GIAO approximation (Cheeseman et al., 1996;Chen et al., 2005). MCIs provide a measure of electron sharing among the atoms considered, whereas delocalization indices (DIs) are two-center electron sharing indices (i.e., two-center MCIs) that measure electron sharing between two atoms (Fradera et al., 1999(Fradera et al., , 2002 and are required to calculate both PDI and FLU. Although several partitions can be used to define the atomic regions needed to calculate DIs and MCIs (Matito et al., 2006), we made use of the molecular partition based on the QTAIM (Bader, 1990(Bader, , 1991. MCI, PDI, and FLU indices were obtained with the ESI-3D program (Matito et al., 2007;Matito, 2014). For MCI, PDI, and HOMA, the larger the values of a given ring, the higher SCHEME 2 | Fragmentation considered in the energy decomposition analysis performed in neutral and dicationic anthracene and phenanthrene. its aromaticity; whereas for FLU, the closer to zero, the more aromatic; and for NICS, the more negative the NICS, the more aromatic the ring. Finally, we have also performed calculations of the anisotropy of the induced current density (ACID) (Herges and Geuenich, 2001;Geuenich et al., 2005) to analyze the ringcurrents at the B3LYP/6-311+G(3df,3pd) level.
The bonding analysis was performed by means of an energy decomposition analysis (EDA), as implemented in ADF (te Velde et al., 2001). Calculations were done at the BLYP/TZ2P level, with the optimized geometries obtained at B3LYP/6-311++G(3df,3pd) level. The bonding analysis was carried out considering two quartet triradicals, fragment 1 (ααα) + fragment 2 (βββ) (see Schemes 2, 3) for the formation of either anthracene or phenanthrene, as well as, of the tetramethylated and dicationic pyrilium pyrano-chromene isomers, whereas for the dicationic anthracene and phenanthrene two triplet diradicals, fragment 1 (αα) + fragment 2 (ββ) were employed. From the defined fragments, the bonding energy was decomposed into two major components (Equation 1; Bickelhaupt and Baerends, 2000;Fernández and Bickelhaupt, 2014;Bickelhaupt and Houk, 2017): In this formula, the strain energy E strain is the amount of energy required to deform the separated triradical or diradical fragments from their equilibrium structure to the geometry that they acquire in the molecule. The interaction energy E int corresponds to the actual energy change when the prepared fragments are combined to form the overall molecule. It is analyzed in the framework of the Kohn-Sham MO model using a Morokuma-type (Kitaura and Morokuma, 1976;Morokuma, 1977) decomposition of the bonding energy into electrostatic interaction, exchange (or Pauli) repulsion, orbital interactions, and dispersion forces (Equation 2; Ziegler and Rauk, 1977): The term V elstat corresponds to the classical electrostatic interaction between the unperturbed charge distributions of the prepared (i.e., deformed) fragments and is usually attractive. The Pauli repulsion, E Pauli , comprises the destabilizing interactions between occupied MOs and it is responsible for the steric repulsion. The orbital interaction, E oi , accounts for bond pair formation (interactions between singly occupied orbitals on one moiety with singly occupied orbitals on the other), charge transfer (donor-acceptor interactions between occupied orbitals of one fragment with unoccupied orbitals of the other, Frontiers in Chemistry | www.frontiersin.org SCHEME 3 | Fragmentation considered in the energy decomposition analysis performed in the pyrano-chromenes derivatives of anthracene and phenanthrene. including the HOMO-LUMO interactions), and polarization (empty-occupied orbital mixing on one fragment due to the presence of another fragment). The orbital interaction energy can be decomposed into the contributions (σ and π contributions in our systems, i.e., E oi = E oi σ + E oi π ) from each irreducible representation Γ of the interacting system using the extended transition state (ETS) scheme developed by Ziegler and Rauk (1977). Finally, the E disp term takes into account the interactions that are due to dispersion forces.
In the EDA, open-shell fragments were treated with the spin-unrestricted formalism but, for technical reasons, spinpolarization was not included. This error causes the studied bond to become in the order of a few kcal·mol −1 too strong. To facilitate a straightforward comparison, the EDA results were scaled to match exactly the regular bond energies (the correction factor applied to all components of the EDA is consistently in the range 0.95-0.99 in all model systems and does therefore not affect trends).

RESULTS AND DISCUSSION
This section is organized as follows. First, we discuss the geometry and isomerization energy of the eight systems analyzed. Second, we quantify the local aromaticity of their rings. Third, we perform energy decomposition analyses using the fragments as depicted in Schemes 2, 3. Finally, we discuss the H···H interaction in the bay region of neutral and cationic phenanthrene and phenanthreneO2 in the framework of the QTAIM theory.
Whereas, single and double bonds are mostly alternated in the three rings of both neutral anthracene and phenanthrene, when two electrons are removed (Figure 1) there is a reorganization, with some bonds elongated and others shortened in each isomer, with a localization of the π system mostly in the terminal rings. In particular, for dicationic anthracene, the π bonds become more localized in the two terminal rings in an alternated but delocalized way, whereas the C-C bonds of the central ring get more single-bond character, at variance with the neutral isomer (see schematic representation in Figure 1). On the other hand, for dicationic phenanthrene, the double bonds are localized in both terminal and central rings, but now the alternation (i.e., delocalization) is only observed for terminal ones, like for charged anthracene.
For tetramethylated-pyrano-chromenes, the methylation in the terminal rings causes a breaking of conjugation (Figure 2), with a unique double bond localized in each terminal ring. Meanwhile, the central ring in both isomers does present both single and double bonds. On the other hand, pyrilium derivatives again present conjugation along the three rings in both conformations (Figure 2).

Aromaticity
The above discussion on the change in geometry and IPs when two electrons are removed from either anthracene or phenanthrene is directly connected to the aromaticity of these systems. Several aromaticity indices have been calculated, i.e., geometrical HOMA, electronic FLU, PDI, and MCI, and magnetic NICS. We will focus our discussion on the MCI results, as it has been previously shown in several tests undertaken that presents the best behavior among different indicators of aromaticity (Feixas et al., 2008(Feixas et al., , 2010. Table 1 encloses the values of all the aromaticity criteria. In the case of NICS, we list only NICS(0) in Table 1 because NICS(0) zz , NICS(1), and NICS(1) zz given in Table S2 provide the same trends. First, in case of neutral isomers, the terminal rings of phenanthrene (MCI = 0.046) are clearly more aromatic than the central one (MCI = 0.018). The five aromaticity criteria point out in the same direction, thus supporting the two Clar π-sextets localized in the two terminal rings (Scheme 1). On the other hand, the aromaticity criteria show a similar aromaticity for terminal and central rings in anthracene (MCI = 0.029 and 0.027, respectively). In particular, HOMA, FLU, and NICS assign a somehow larger aromaticity to the central ring, whereas the inner and outer rings have an almost equivalent aromatic character according to both electronic-based criteria (PDI and MCI). This result agrees with its migrating Clar's π-sextet (Scheme 1). By comparison, terminal rings are more aromatic in phenanthrene, whereas the central one is more aromatic in anthracene (Szczepanik et al., 2018).
Considering the Platt's perimeter model (Platt, 1949), neutral isomers follow 4n+2 Hückel's rule, with 14 π electrons, whereas the charged ones does not, because they have 12 π electrons, thus they should be considered antiaromatic based on this model. Then, once two electrons are removed, dicationic phenanthrene shows a clear reduction of aromaticity in both types of rings (MCI = 0.019 and 0.010 for terminal and central rings, respectively). Despite their reduction, terminal rings continue being somewhat aromatic, as proven by all criteria, except of NICS. The decrease of aromaticity also happens for the central ring, which might be considered non-aromatic. On the other hand, dicationic anthracene presents a different behavior: the terminal (MCI = 0.028) rings are now equally (PDI and MCI) or slightly more aromatic (HOMA and FLU) than in the neutral species, whereas the aromaticity of the central ring (MCI = 0.013) is slightly reduced. Both terminal and central rings of dicationic anthracene are more aromatic than for dicationic phenanthrene. Further, in both cases, the terminal rings are more aromatic than the central one, also in agreement with what is expected based on Clar's rule (Scheme 1). ACID plots (Figure S1) support the aforementioned changes in aromaticity when going from neutral to cationic phenanthrene and anthracene. Anthracene shows diamagnetic ring-currents in the inner and outer rings that become somewhat more intense in the outer rings when two electrons are removed. On the other hand, the ACID plot of phenanthrene indicates more aromatic outer than inner rings that lose aromaticity in the dicationic form. However, it must be pointed out the dependence of ACID on the location of the integration plane (Fliegl et al., 2016).
However, at variance with the prediction of the Clar model of Scheme 1, Voronoi deformation density (VDD) charges in Figure 3 show that, when the two electrons are removed from both anthracene and phenanthrene, the atoms of the terminal rings experience the largest loss of electronic charge (more positively charged atoms). More importantly, the Clar valence structure of anthracene dication corresponds to two benzyl cations located in para position, whereas that for phenanthrene dication corresponds to a vicinal dication ( Figure 3B). These expected models are supported by the VDD charges ( Figure 3A) with +87 me − for carbons in para in anthracene +2 and +71 me − for carbons in ortho in phenanthrene +2 . The positive charges in ortho cause an important electrostatic repulsion that is reduced by partial delocalization of the positive charge to outer rings, thus decreasing the aromaticity of these rings. In this way, we can justify the loss of aromaticity of the outer rings of phenanthrene when two electrons are removed.
From the aromaticity analysis, in case of neutral systems the higher stability of the kinked isomer correlates with its larger aromaticity. On the other hand, for charged systems, despite both isomers present a closer aromatic character; anthracene dication is more aromatic, which also correlates with its larger stability. These results are in agreement with previous work by Dominikowska and Palusiak (2011). Furthermore, these trends in aromaticity also explain the ionization potentials discussed in the previous section, being smaller for anthracene as the loss of two electrons in this case transforms the less aromatic neutral system into the most aromatic cationic benzenoid.
With respect to the third set of isomers, i.e., tetramethylatedpyrano-chromene derivatives, the dimethylation of each terminal ring causes these rings to be non-aromatic (MCI = 0.002 for both isomers), whereas the central ring is clearly aromatic (MCI = 0.039 and 0.042 for anthraceneO2 and phenanthreneO2, respectively). All aromaticity criteria point out a similar aromaticity for both isomers, which does not explain the larger stability of the linear one. So, further analysis is required in order to understand the change in the relative stability of the linear and kinked isomers when we move from neutral to dicationic species. Finally, the behavior is opposite for the pyrilium derived systems, as in this case phenanthrenePyr is slightly more     stable than anthracenePyr, which is likely to be because the differences in aromaticity are lower for these systems than for the couple anthracene/phenanthrene. In particular, the outer rings in phenanthrenePyr with MCI = 0.020 are significantly less aromatic than their analog rings in phenanthrene (MCI = 0.046).

Energy Decomposition Analysis
The larger stability of phenanthrene than anthracene, and the opposite for their corresponding dicationic systems, can be understood by means of a quantitative Kohn-Sham MO analysis complemented with an energy decomposition analysis (EDA, Table 2). We have to mention that the EDA of the anthracene and phenanthrene was already discussed in detail in our previous work (Poater et al., 2007b). Here, we briefly summarize the main results and provide the values obtained at the BLYP/TZ2P//B3LYP/6-311++G(3df,3pd) level of theory. Both phenanthrene and anthracene can be built from two identical 2-methtriyl-phenyl fragments (Scheme 2) in their quadruplet state. The three unpaired electrons are located in three singly occupied molecular orbitals (Figure 4). These are the σ A and σ S orbitals in the σ electron system and the π orbital in the π electron system. The two former are antisymmetric and symmetric regarding the sign of their largeamplitude lobes. The more favorable total bonding energy for phenanthrene ( E = −4.7 kcal mol −1 ) is due to the interaction energy ( E int = −4.8 kcal mol −1 ), as the deformation of the fragments is almost the same in both cases ( E strain = 0.1 kcal mol −1 ). E int can be further decomposed as in Equation (2). The more favorable E int is mainly attributed to Pauli repulsion ( E Pauli = −16.1 kcal mol −1 ), whereas both attractive electrostatic and orbital interaction terms are more stable in anthracene ( V elstat = 9.0 and E oi = 2.9 kcal mol −1 ). However, it has been previously proven that this observation is misleading, as it does not only reflect the pure changes in bonding because of "flipping" the two fragments to form either anthracene or phenanthrene, but it also contains the effect of structural relaxation that is induced by the changes in intrinsic bonding which also masks the latter (Bickelhaupt et al., 2002;Bickelhaupt and Baerends, 2003;Poater et al., 2007b). It has been shown that a system with larger Pauli repulsion may lead to an equilibrium structure with a longer bond in which E Pauli is lower than in the equilibrium structure of the sterically less demanding equivalent system. Thus, the differences between anthracene and phenanthrene cannot be only analyzed based on their equilibrium geometries, as we may get erroneous conclusions. To overcome this issue, the bonding analysis in a deformed phenanthrene has been performed by taking the same fragments as in anthracene, but just flipped, i.e., turn-upside-down approach (El-Hamdi et al., 2011), with the two formed C-C bonds kept fixed to the equilibrium distance of anthracene (formed C-C bonds: 1.399 Å). This modified phenanthrene (phenanthrene Ant Ant ) is still more stable than anthracene ( E = −1.2 kcal mol −1 ), although less than phenanthrene in its equilibrium geometry. The only difference between anthracene Ant Ant and phenanthrene Ant Ant is the connection of the two fragments, i.e., the topology. Phenanthrene Ant Ant presents almost the same steric repulsion between the two fragments as anthracene Ant Ant ( E Pauli = −0.6 kcal mol −1 ), so it is not the determinant term any longer. The difference is mainly due to the π orbital interaction ( E oi π = −1.5 kcal mol −1 ). As can be seen in Figure 4, the π SOMO of the 2-methtriyl-phenyl fragment has small coefficients in the ring and a large one on the exocyclic carbon. A better π-π overlap is achieved when the two exocyclic lobes overlap and this occurs in phenanthrene. On the other hand, E oi σ = 0.9 kcal mol −1 and this result is important because if the stabilizing H···H bonding in the bay region of phenanthrene would exist, as previously stated based on QTAIM, then E oi σ would favor phenanthrene topology.
The next step consisted of using the same two fragments, but now with the formed C-C equal to those of phenanthrene (formed C-C bonds: 1.357 and 1.457 Å). This new geometry (phenanthrene Phe Ant ) already gives a lower Pauli repulsion ( E Pauli = −8.8 kcal mol −1 ), but more importantly, together with a larger stabilization of the π orbital interaction ( E oi π = −2.6 kcal mol −1 ). The final step involves going to the fully relaxed phenanthrene, which clearly shows an even increased reduction of E Pauli , as stated above, due to the bending of the C-H bonds in the bay region. So, based on these model systems, the more efficient bonding in the π-electron system in phenanthrene is responsible of its larger stability compared to its linear isomer, which is supported by its larger aromaticity. Such larger E oi π for phenanthrene is also supported by the better overlap between the π SOMO of each fragment (0.243 vs. 0.196 for relaxed phenanthrene and anthracene, respectively; in case of phenanthrene Ant Phe and phenanthrene Ant Ant , the overlaps are 0.236 and 0.230, respectively). Now, as stated above, dicationic anthracene becomes more aromatic than charged phenanthrene, especially due to the decrease of aromaticity of the latter. However, despite their closer aromatic character, dicationic anthracene is more stable than phenanthrene by 16.2 kcal mol −1 ( E), with a larger difference than between the neutral isomers (see Table 3). From the EDA analysis, the strain energy is smaller for phenanthrene dication, which gives a E int = 22.2 kcal mol −1 in favor of dicationic anthracene. Like in the previous case, Pauli repulsion is smaller in TABLE 3 | Energy decomposition analysis (EDA) data for dicationic anthracene and phenanthrene (in kcal mol −1 ), and the corresponding turn-upside-down derivatives. charged phenanthrene ( E Pauli = −64.7 kcal mol −1 ), however now it cannot compensate the more stable electrostatic and orbital interaction terms ( V elstat = 31.5 and E oi = 55.8 kcal mol −1 ) of charged anthracene. Electrostatic interactions for the charged species are smaller than for the neutral anthracene and phenanthrene, as expected from the fact that the two approaching fragments in charged systems bear a positive charge. Interestingly, in dicationic anthracene and phenanthrene both the E oi σ and E oi π favor the straight isomer. As done before, the turn-upside-down approach has also been applied for these species to better understand the larger stability of dicationic anthracene vs. its kinked isomer. A deformed phenanthrene from the same fragments as for anthracene is built, first with the same new C-C bonds as those of dicationic anthracene (formed C-C bonds: 1.413 Å, phenanthrene 2+ Ant Ant ), and next with the new C-C bonds as those of relaxed dicationic phenanthrene (formed C-C bonds: 1.445 and 1.475 Å, phenanthrene 2+ Ant Phe ). Again, phenanthrene 2+ Ant Ant and anthracene 2+ Ant Ant differ only in the topology. Phenanthrene 2+ Ant Ant is less stable than relaxed anthracene by 24.5 kcal mol −1 . Pauli repulsion slightly favors the anthracene dication ( E Pauli = 4.5 kcal mol −1 ), but the main contribution comes from more stabilizing orbital interactions ( E oi = 18.3 kcal mol −1 ) both σ-( E oi σ = 0.9 kcal mol −1 ) and π-orbital interactions ( E oi π = 17.4 kcal mol −1 ). Because the SOMOs σ are more or less the same in the neutral and charged species (Figure S2), the E oi σ stabilization is basically equal (difference of only 0.2 kcal mol −1 ). For the π-system, the situation is totally different (Figure S3). In the neutral systems, there is a SOMO π that interacts better in phenanthrene topology than in anthracene one favoring the kinked isomer. In charged species, the fragments contain 6π electrons in three occupied π-MOs. Now, the E oi π contains neither bond pair formation nor charge transfer interaction, only polarization. We consider that polarization is more favorable in the linear topology because positive charges are located mainly in the outer rings and these rings are more separated in the anthracene dication ( Figure S4). Once the new C-C distances are changed to those of relaxed charged phenanthrene (phenanthrene 2+ Ant Phe ), the terms are already very close to those of fully relaxed dicationic phenanthrene. So, despite Pauli repulsion is reduced, electrostatic and orbital interaction terms are less attractive in phenanthrene 2+ Ant Phe than phenanthrene 2+ Ant Ant . Overall, our results prove that also for dicationic isomers, the E oi π term plays a critical role in determining the relative isomeric energies. This term favors the kinked structure for neutral and the linear for dicationic anthracene and phenanthrene.
Next, we analyze the pair of kinked and straight tetramethylated-pyrano-chromene isomers having both only one single Clar sextet. These isomers can also be formed from the same fragments (Scheme 3). AnthraceneO2 presents a lower total binding energy than for phenanthrene isomer by 6.2 kcal mol −1 ( E), from which 3.5 kcal mol −1 comes from the strain energy and 2.7 kcal mol −1 from the interaction energy. Further, again the phenanthrene derivative has a lower Pauli repulsion ( E Pauli = −33.8 kcal mol −1 ), which is compensated by more attractive electrostatic and orbital interactions in anthraceneO2 ( V elstat = 18.9 and E oi = 18.6 kcal mol −1 ). In this case, the more stable E oi for anthraceneO2 mainly comes from the σ instead of the π interactions between the MOs of the two fragments, but like in the dications both favor the linear isomer. The E oi π = 4.2 kcal mol −1 favoring the anthraceneO2 isomer, proves its determinant role at assigning the larger stability to either the kinked acene (neutral) or the linear one (dicationic and tetramethylated-pyrano-chromene derivatives). As for the neutral systems, this more favorable E oi π for anthraceneO2 is supported by the overlaps between the π SOMOs of each fragment (0.228 vs. 0.247 for relaxed phenanthreneO2 and anthraceneO2, respectively). For phenanthreneO2 and anthraceneO2, the overlaps between the π SOMOs do not differ much in the two topologies because in the AntO fragment the contribution of the π SOMO in the ring and in the exocyclic carbon are very similar (see Figure S5).
Likewise, we applied the turn-upside-down approach to tetramethylated-pyrano-chromene derivatives. Unfortunately, it cannot be applied straight because the presence of the oxygen in the terminal rings causes that the fragment does not have the C2 axis to get the correct fragment when turned-upside-down in order to build the other isomer. So, TABLE 4 | Energy decomposition analysis (EDA) data for tetramethylated-pyrano-chromenes (in kcal mol −1 ), and the corresponding turn-upside-down derivative. instead of breaking anthraceneO2 system in two fragments (AntO + AntO') and rotating one of them, we need twice the same fragment (AntO), but one rotated, to build phenanthreneO2 AntO AntO isomer (formed C-C bonds: 1.392 Å) and phenanthreneO2 AntO PheO (formed C-C bonds: 1.393 and 1.413 Å) isomer, as above (Table 4). PhenanthreneO2 AntO AntO is 7.0 kcal mol −1 less stable than anthraceneO2, which again is due to less stabilizing electrostatic and orbital interactions, which cannot compensate its lower Pauli repulsion. Then, in case of phenanthreneO2 AntO PheO , the values are already very close to phenanthreneO2, and the behavior very similar to that of dicationic systems, i.e., despite Pauli repulsion is reduced, electrostatic and orbital interaction terms are less attractive. Just for comparison, the overlaps between the π SOMOs in case of phenanthreneO2 AntO PheO and phenanthreneO2 AntO AntO are 0.226 and 0.230, respectively, to be compared to the discussed above of 0.247 for anthraceneO2.
For the last couple of systems, i.e., anthracenePyr vs. phenanthrenePyr (Table 5 and Scheme 3), the latter appears to be only 0.9 kcal mol −1 than the former ( E). This difference arises due to more favorable E int term of the kinded system ( E int = −3.2 kcal mol −1 ), despite its unfavorable strain energy ( E strain = 2.2 kcal mol −1 ). At difference with the previous couple of pyrano-chromene systems, now the more favorable V elstat and E oi of the linear system cannot compensate its much larger E Pauli , and, for this reason, phenanthrenePyr is more stable. From the turn-upsidedown approach we can observe that when phenanthrenePyr is built from the same fragments of anthracenePyr, i.e., phenanthrenePyr AntPyr AntPyr or phenanthreneO2 AntPyr PhePyr , Pauli repulsion increases in phenanthrenePyr and causes the kinked system to be less stable, proving the determinant role of E Pauli . At variance with what is observed in phenanthrene, the E oi π term in phenanthrenePyr AntPyr AntPyr is more favorable for the linear than for the kinked isomer. Therefore, π-interactions are stronger in the most stable linear topology. However, relaxation of phenanthrenePyr AntPyr AntPyr to phenanthrenePyr generates a slightly more stable kinked structure because of the large reduction in Pauli repulsions.

QTAIM
The QTAIM atomic energies of the H atoms involved in the H···H interaction in the bay region of phenanthrene were calculated for the neutral and dicationic systems (see Figure 5). In case of the neutral isomers, the H atoms in the bay region of phenanthrene are stabilized by −5.2 kcal mol −1 compared to the non-interacting hydrogen atoms in linear anthracene. This stabilization was used to justify the larger stability of the kinked isomers (Matta et al., 2003), as stated above. However, in case of the dicationic systems, the same H atoms in phenanthrene dication are even more stabilized, by −11.2 kcal mol −1 , compared the non-interacting hydrogen atoms in dicationic anthracene. This should translate into a more stable charged phenanthrene than anthracene, which is not the case, thus providing an additional evidence that QTAIM atomic energies of few selected atoms cannot be used to justify isomerization energies. In addition, the same trend is found for the tetramethylated-pyrano-chromene derivatives, for which the bay H atoms in the kinked isomer are −3.0 kcal mol −1 more stable than in the linear one. Again, this result does not concur with the greater stability of the linear isomer.

CONCLUSIONS
The main conclusion of this work is that the relative stability of kinked vs. straight topologies in polycyclic aromatic hydrocarbons (PAHs) depends on the aromaticity of the rings of the isomers. So, the answer to the question whether kinked is more stable than straight depends on the aromaticity of the two isomers. Phenanthrene is more stable than anthracene due to the larger stability of the π-system of the former, which is more aromatic. When two electrons are removed, i.e., dicationic systems are analyzed, the reverse trend is obtained, so the linear isomer is more stable than the kinked one. The larger stability of the former also correlates with its larger aromaticity, despite the difference in aromaticity between the isomers is now smaller. Finally, tetramethylated-pyrano-chromene isomers almost present the same aromatic character, whereas pyrilium derived systems show divergent aromaticities between central and terminal rings, i.e., central is more aromatic for the linear and terminal for the kinked. In these cases, when no isomer is clearly more aromatic, the relative stability is determined by the balance between the lower steric repulsion of the kinked vs. the more favorable electrostatic and orbital interactions of the linear.
Just to conclude, our study together with other works (Poater et al., 2006a(Poater et al., ,b, 2007a(Poater et al., ,b, 2017 clearly contradicts previous studies based on QTAIM analysis stating that the larger stability of phenanthrene is due to the stability of the H···H bridging interaction (Matta et al., 2003). This interaction is more stable for the kinked isomer in the three sets of isomers analyzed; however only in one case, the kinked system is more stable.

AUTHOR CONTRIBUTIONS
JP and MS conceived and designed research. JP, MD, and MS analyzed data, interpreted results and wrote the paper. JP performed the calculations.

ACKNOWLEDGMENTS
We are very grateful for the advice given by Dr. Ouissam El Bakouri (Uppsala University) on the calculations of the ACID graphics. This work has been supported by the Ministerio de Economía y Competitividad (MINECO) of Spain (Projects CTQ2017-85341-P, CTQ2016-77558-R, and MDM-2017-0767) and the Generalitat de Catalunya (project 2017SGR39, Xarxa de Referència en Química Teòrica i Computacional, and ICREA Academia 2014 prize for MS). The EU under the FEDER grant UNGI10-4E-801 (European Fund for Regional Development) has also funded this research.