Molecular asphaltene models based on Clar sextet theory

Abstract

Asphaltenes, the base components in asphalt used for pavements and roads, are high-molecular-weight polycyclic aromatic hydrocarbons with long aliphatic chains and polar functionalities. Significant efforts have been made to develop representative asphalt models for molecular simulations, and several molecular structures have been proposed for asphaltenes based on experimental results and theoretical predictions. Many properties, such as viscosity, depend on the underlying detailed chemical structure of asphaltenes, since it defines the interactions at the nanoscale that dominate the formation of nanoaggregates. Among these interactions, π–π stacking plays a central role in the formation of clusters at the nanoscale level, with strong dependence on the electronic arrangement of π electrons in the polycyclic aromatic core of asphaltenes. Here, based on Clar sextet theory and density functional theory (DFT) calculations, we propose asphaltene molecules that optimize the distribution of π electrons in the polycyclic aromatic core and minimize geometrical strain within the molecule. Accordingly, more stable isomers than the asphaltene molecules used in current asphalt models are proposed. These new isomers include a more representative arrangement of the π electrons, crucial for further description of intermolecular interactions and formation of nanoaggregates. Considering these basic guidelines in future molecular modeling of asphalt will improve the development of more accurate models that describe the behavior of these complex materials at the nanoscale.

---

1. Introduction

Although most of the properties of asphalt materials are known to depend on their underlying chemical and nanoscale structure,¹ the detailed molecular composition of asphalts remains unknown in many cases. The uncertainty arises because a large fraction of crude oil distillation includes millions of distinct organic molecules.¹ This complexity makes it extremely difficult to elucidate the precise chemical structure for all its constituents. Solubility differences and chromatography studies make it possible to categorize the constituents into four different classes, commonly described as saturates, aromatics (or naphthene aromatics), resins (also known as polar aromatics) and asphaltenes, which are usually referred to as SARA.² Beyond classification, the SARA model also provides an overall description of the chemical functionalities that might be present in these molecules. Many intermolecular interactions depend on the functional groups and also on their specific location within the molecule, e.g. the arrangement of aromatic rings in a polycyclic aromatic system. Therefore, describing the chemical structure accurately is essential to reproduce nanoscale and bulk properties of asphalt materials.

Among these different types of molecules in the SARA model, asphaltenes deserve special attention.¹ These are high-molecular-weight molecules, insoluble in light alkanes but soluble in benzene and toluene, whose structure and interactions at the nanoscale strongly affect rheological and mechanical properties of asphalt.^3,4 Its fused aromatic ring (FAR) region is similar to a polycyclic aromatic hydrocarbon (PAH), but long aliphatic chains, and some polar functionalities are also present in asphaltenes. Being similar to PAHs, the aromaticity that results from the delocalization of π electrons in the FAR region is responsible of π–π interactions between asphaltene molecules, which lead to the formation of clusters at the nanoscale.³ Therefore, understanding the implications of different geometries and sizes on the stability of the molecules is essential to understand their supra-molecular arrangement.

The central role of asphaltenes in the use of petroleum resources, together with the uncertainty about their chemical structure brings an interesting challenge for molecular modeling and simulation. In this regard, considerable effort has been devoted to elucidating the detailed chemical structure and the supra-molecular arrangement of asphaltenes. Numerous models have been put forward, which can be categorized into island^3,5 and archipelago^6–8 models, both supported by experimental studies.^9,10 Island models support the idea of asphaltene molecules composed by one medium-size FAR region with aliphatic chains. Archipelago models propose asphaltene molecules where many FAR regions are connected by short aliphatic carbon chains. Debate arises also in regards to the supra-molecular arrangement, and both colloidal distribution of molecular aggregates and molecularly disperse models have been reported.^2,5,11 Concerning the development of asphalt models, a large body of remarkable work has been done by Greenfied,^4,12–14 Mullins,^3,5,15 and others.^16,17 Mullins modified the previously-proposed Yen model,¹⁸ by suggesting new asphaltene chemical structures based on the most probable molecular weight within island architecture, and also describing the formation of nanoaggregates and clusters.^3,5 More recently, Li and Greenfield proved that some of the molecules currently used in these asphalt models present high internal energies¹⁹ due to the so-called “pentane effect”.^19,20 This effect occurs in asphaltene molecules due to the steric repulsion between the carbon aliphatic chains, causing the FAR region to emerge out of the plane. The geometry deformation caused by this effect increases the energy of the molecule. Based on density functional theory (DFT) calculations, Li and Greenfield proposed new isomers for the molecules included in the Yen–Mullins model, which do not undergo pentane effect and constitute a further improved model for asphaltenes.⁴ These new isomers are included in a new complete model of asphalt that should serve as reference for future work in the field.⁴

While the latter work place a significant emphasis on geometrical strain within the molecule to avoid unphysical geometries and high-energy isomers, it does not optimize the electronic arrangement in the FAR region of the molecule. It should be noted that due to the importance of π–π stacking in the formation of clusters at the nanoscale level, optimizing the distribution of π electrons in the FAR region is crucial for further description of the interactions and properties. Therefore, this paper takes a further step beyond internal stress when proposing new asphaltene isomers with a more representative arrangement of aromatic rings in the system.

Furthermore, this paper outlines the efficacy of Clar's theory of the aromatic sextet^21–23 for predicting the stability of aromatic systems. The application of this theory to asphaltene FAR regions has been performed thoroughly by Ruiz-Morales et al.,^1,9,24–26 fitting its predictions to experimental data, and suggesting 4 to 10 fused aromatic rings. This framework allows defining very simple but elegant rules towards proposing asphaltene molecules with optimized distribution of π electrons in the FAR region (more detailed description on Clar sextet theory is provided in ESI†). Accordingly, using this theory, isomers that are more stable than those in current asphalt models are proposed. If one considers that asphaltenes are essentially stable for geological time, less stable aromatic structures (as determined by DFT calculations) are those structures which have lower occurrence probability,¹⁹ and therefore, lower-energy chemical structures provide more realistic molecules for molecular simulations. Therefore, we aim to suggest alternative molecular structures for asphaltenes by incorporating Clar sextet theory rules to the study of current asphaltene models, complementing earlier work. The effort to elucidate the relationship between mechanical properties and chemical structure, named as chemo-mechanics,⁴ will be only possible if we build this understanding in a robust bottom-up approach. Therefore, in this paper, standing on Clar's theory of the aromatic sextet and DFT calculations, we propose more stable chemical structures for some of the asphaltene molecules used in present asphalt models in order to improve their accuracy. Beyond that, we suggest some general guidelines to be considered for the molecular modeling process. This strategy helps to propose asphaltene isomers with optimal aromaticity distribution in the polycyclic aromatic core and it is applicable to any asphaltene molecule belonging to either island or archipelago models.

2. Methods

2.1 Computational details

DFT calculations are performed with Gaussian 09 package.²⁷ Molecular geometries are fully optimized (see ESI† for XYZ coordinates of optimized geometries) using DFT with Becke's three-parameter hybrid exchange functional and Lee–Yang–Parr correlation functional, commonly referred to as B3LYP.^28,29 This is a hybrid functional that goes beyond the generalized gradient approximation (GGA) by including some exact exchange energy calculated from Hartree–Fock theory. Additionally, 6-31+G* basis set is used for both the optimization and electronic structure calculations. This presents a split-valence double-zeta basis set with both polarization and diffuse functions included. Furthermore, frequency calculations have been performed on the optimized geometries. No negative frequencies have been found and therefore it is proved that we are dealing with ground state structures (see ESI† for details concerning frequencies).

For comparing the planarity of the FAR region of the different molecules, the rings are colored attending to Cremer–Pople pucker amplitude values, which is a puckering spherical coordinate system initially described for six-membered pyranose rings,³⁰ but extensible to any other ring with three or more members. In this coordinate system the radius Q means the magnitude of puckering, measuring the deviation from the perfectly flat six-membered ring (Q = 0). The values of Q obtained by this method are color-coded^31,32 by the PaperChain representation implemented in VMD software.³³ More information about this methodology is included in ESI.†

2.2 Clar's theory of the aromatic sextet

Clar's theory of the aromatic sextet^21–23 is used for analyzing the FAR region of asphaltene molecules. This theory states that the most important representation of a PAH is one having the maximum number of π-sextets, depicted by inscribed circles, and minimum number of fixed double bonds. PAHs with more sextets are more stable due to the delocalization of π electrons within these resonant sextets. Therefore, different isomers can have different number of resonant π-sextets due to the different geometrical arrangement of the fused rings. Isomers with larger number of resonant sextets would be more stable. Based in these simple but useful ideas, we proposed new isomers for some asphaltene molecules used in current models of asphalt. More information about Clar's theory of the aromatic sextet is provided in ESI.†

3. Results and discussion

DFT calculations were performed in a small set of asphaltene molecules commonly employed in current asphalt models. The structure of these molecules is analyzed using Clar sextet theory in order to propose lower-energy isomers. The isomers proposed by these means present an optimized arrangement of π sextets of electrons in the FAR region. Fig. 1 shows the chemical structure of the asphaltene molecules under consideration, which we refer to as A1, A2 and A3 throughout the following discussion. These three molecules have been selected for being the most commonly employed structures in current asphalt models. A1 has been the asphaltene molecule employed by Greenfield et al. in many of their previous models of asphalt.^12,13,34,35 This molecule presents a non-planar FAR core due to the existence of pentane effect between two rings in the structure. This effect has been recently described for asphaltenes,¹⁹ although the specific case of A1 molecule has not been explored yet. Since A1 has been used in many models, we find it reasonable to discuss the geometrical issues of this structure as well, together with the improvement associated to its π-sextet distribution. It has to be mentioned that the kind of steric repulsion discussed in this case does not involve any aliphatic chains, as in the so-called pentane effect and mainly focuses on steric interactions between aromatic rings in the FAR core. Additionally, we consider A2 and A3 molecules, which are some of the latest asphaltene molecules proposed also by Li and Greenfield when improving Yen–Mullins' model.^4,19 Studying these two molecules allows us to demonstrate the importance of both geometrical issues and π-sextet distribution when studying asphaltene chemical structures. While A2 presents a chemical conformation that is also susceptible to steric repulsion between aromatic rings, A3 does not present any geometrical issues. However, it can still be improved following Clar sextet theory.

Fig. 1 Asphaltene molecular structures considered in the present study. (a) Asphaltene from previous Greenfield models^12,13,34,35 and (b), (c) Greenfield-modified asphaltenes from the Yen–Mullins model.⁴

In the first step, aiming to reduce the size of the model systems for DFT calculations, we consider smaller molecules that contain the same FAR core, which means reducing asphaltene molecules to the PAH that represents them. These PAH molecules would retain the same electronic distribution and geometrical effects in the aromatic system, but with considerably smaller number of atoms. To confirm this assumption, three different molecular structures for representing A1 are investigated using Cremer–Pople pucker amplitude to evaluate ring deformation. A complete description of these three systems is provided in ESI.† The results support that including the aliphatic chains in the calculations is not needed for the purpose of this study, since the deformation of the different rings is almost the same in the three models. Therefore, a first step in the designing process of asphaltene molecules is reducing the system to a PAH where no aliphatic chains are included. Doing so, the computational cost is reduced and the design will focus on improving the geometrical arrangement of aromatic rings in the system. Accordingly, in the following discussion we use these simplified molecules to represent asphaltenes. We refer to them as A1, A2 and A3. If at some point the complete molecule, with aliphatic chains, is discussed we will specify this.

Using these simplified structures, we study the arrangement of aromatic rings in the system to propose new more-stable isomers. Fig. 2 shows the rational design of this new arrangement of aromatic rings in the A1 molecule. A1-M1 is the original molecule, while A1-M2 and A1-M3 stand for the two new proposed isomers. Fig. 3 shows the side view of the optimized molecular geometries together with the energies relative to the most stable isomer, calculated at B3LYP/6-31+G* DFT level of theory. The top part of Fig. 2 shows the chemical structure of the three molecules under consideration, where aromatic rings undergoing pentane effect are highlighted in red. In the A1-M3, the ring highlighted in blue indicates that this is the ring that has changed its position from the previous isomer A1-M2. Besides, the bottom part of Fig. 2 shows the Clar configurations for these three isomers, where cycloalkane rings that do not take part in the Clar formula are shadowed in grey for clarity. The main purpose of this step of the design process is to rearrange the FAR region in such a way that the number of Clar π-sextets is maximized, while avoiding any possible pentane effect in the geometry. As it can be noticed, A1-M1's original molecule already undergoes pentane effect between the two rings in red. This steric interaction leads the geometry of the FAR core to emerge out of the plane, as shown in Fig. 3a. In A1-M2, considering a different arrangement of aromatic rings where the most external ring involved in the steric repulsion in A1-M1 is now placed in a position that avoids such interaction, leads to a 99.58 kJ mol⁻¹ less energy. At the same time, A1-M2 increases the number of π-sextets in the structure. Actually, the Clar configuration of A1-M1 (Fig. 2) leads us to the possible positions where an additional π-sextet could be located, guiding the design of the new arrangement of rings in A1-M2. Thus, by following the rules of Clar sextet Theory, A1-M2 isomer increases the number of Clar sextets (from 3 Clar rings to 4) lowering the energy; the latter arrangement still undergoes certain strain due to the proximity of the rings highlighted in red in the A1-M2 (Fig. 2). This strain is completely reduced in A1-M3 structure, where the new position of the aromatic ring in blue does not interact with adjacent rings. This additional change in the FAR region implies some extra stabilization of 12.22 kJ mol⁻¹. The smaller amount of energy change between A1-M2 to A1-M3 is due to fact that in this case only geometrical issues are involved, and no additional Clar π-sextet is included in the system. Overall, A1-M3 not only does avoid steric repulsion, but also increases the number of aromatic π-sextets in the structure with respect to A1-M1. Accordingly, the higher number of Clar sextets in A1-M2 and A1-M3 provides higher stability to these systems. The total energy difference between the existing model (A1-M1) and the new proposed isomer (A1-M3) is 111.80 kJ mol⁻¹ indicating that A1-M3 isomer is much more stable than existing model. In addition, the changes in energy from A1-M2 and A1-M3 corroborate that the FAR region of A1-M2 still undergoes some geometrical strain due to the proximity of those rings highlighted in red in Fig. 2. Furthermore, Fig. 3 shows that the structure of the molecules becomes more planar as the pentane effect is avoided. This extra planarity comes from avoiding steric repulsion while increasing the aromaticity of the FAR region due to the additional π-sextet of electrons.

Fig. 2 Simplified model of the three different isomers of A1 asphaltene; rings involved in pentane effect are highlighted in red, while relocated rings avoiding such effect are in blue. Bottom part shows the Clar formula of the different isomers, where cyclohexane rings are shadowed in grey to indicate that they do not take part is the Clar formula.

Fig. 3 Side view of the optimized structures of A1, A2 and A3 asphaltene isomers at B3LYP/6-31+G* level of theory. The energy relative to the most stable one is shown in parenthesis (in kJ mol⁻¹). (a) Side view of the optimized structures of A1 asphaltene isomers. (b) Side view of the optimized structures of A2 asphaltene isomers. (c) Side view of the optimized structures of A3 asphaltene isomers.

A similar analysis is followed in the study of A2 asphaltene. As we mentioned before, this molecule was firstly suggested in the Yen–Mullins model,³ and it was later modified by Greenfield et al.⁴ Nevertheless, we still find its molecular structure suitable for further improvement from both geometrical and π-sextet distribution points of view. Again, a simplified PAH model where only hydrogen saturates terminal carbon position is used for the calculations. Fig. 4 shows the new isomers, A2-M2 and A2-M3, proposed for A2-M1. In the top part of the Fig. 4, aromatic rings that undergo pentane effect are highlighted in red. The ring highlighted in blue shows the final position of the ring that has changed its position in the rearrangement process of the FAR region. The bottom part of the figure shows the Clar configurations for these three isomers, where again cycloalkane rings are shadowed in grey for clarity. As shown in the figure A2-M1 isomer undergoes pentane effect, which leads the geometry of the aromatic core to come out of plane (see Fig. 3b) leading to an increase in its internal energy by 62.02 kJ mol⁻¹ compared to A2-M2 isomer and by 65.20 kJ mol⁻¹ compared to A2-M3 isomer; accordingly, the latter in found to be the most stable one among the three isomers. Therefore, moving one of these interacting rings (highlighted in red in A2-M1 isomer) to a non-interacting position (the one attached to the pentagon in A2-M2) decreases the internal energy of the molecule by 62.02 kJ mol⁻¹. This rearrangement of the rings in the PAH is done following Clar rules in such a way that an additional π-sextet is accommodated. However, although the A2-M2 isomer is 62.02 kJ mol⁻¹ lower in energy than A2-M1 due to the presence of an extra π-sextet, it still undergoes some internal steric repulsion due to the interaction of the two rings highlighted in red in the central part of Fig. 4. This kind of interaction is avoided in A2-M3, where the final location of the ring (in blue) does not undergo pentane effect leading to another 3.18 kJ mol⁻¹ energy reduction. Therefore, A2-M2 and A2-M3 isomers increase the number of aromatic π-sextets in the structure from 4 to 5, which is reflected in higher stability, while avoiding steric interactions. Fig. 3b shows the side view of the optimized geometries for the A2 isomers, together with the relative energy in parenthesis.

Fig. 4 Simplified model of the three different isomers of A2 asphaltene; rings involved in pentane effect are highlighted in red, while relocated rings avoiding such effect are in blue. Bottom part shows the Clar formula of the different isomers, where cyclo-hexane rings are shadowed in grey to indicate that they do not take part is the Clar formula.

Besides, we examine A3 asphaltene. This molecule was also suggested in the Yen–Mullins model,³ and later modified by Greenfield et al.⁴ (Fig. 5). A3 is especially interesting from the point of view of the Clar sextet analysis, since it does not undergo pentane effect and therefore the improvement should be exclusively based on the location of additional Clar π-sextets in the structure. The simplified PAH model is used also for the calculations. Fig. 5 shows both the original A3-M1 molecule and a new isomer A3-M2, whose chemical structure is based on the guidelines suggested by Clar's theory of the aromatic sextet. The ring highlighted in red in the top part of the figure is the one consider for relocation, while the blue ring indicates the final location of this ring. The bottom part shows the Clar configurations for both molecules, where cycloalkane rings are shadowed in grey. Moving a single aromatic ring allows the location of an extra π-sextet in the structure, which provides additional stability to the system. Indeed, the optimized geometries at B3LYP/6-31+G* level of theory show that this extra aromatic sextet implies almost 26 kJ mol⁻¹ less internal energy. It is remarkable that in this case there is no pentane effect, and the additional stability is due exclusively to the presence of an extra π-sextet in the system, which enhances the importance of the Clar's π-sextet rules on proposing new asphalt isomers. Both A3-M1 and A3-M2 are planar (see Fig. 3c) and there is no contribution to the stabilization resulted from any geometry deformation, but only from the electronic arrangement in the FAR region.

Fig. 5 Simplified model of the two different isomers of A3 asphaltene; Ring that is suitable of relocation according to Clar's sextet rule is highlighted in red, while relocated ring is in blue. Bottom part shows the Clar formula of the isomers, where cyclo-hexane rings are shadowed in grey to indicate that they do not take part is the Clar formula.

Table 1 summarizes the relative energies discussed throughout the paper, all of them calculated with B3LYP/6-31+G* level of theory. The table shows how A1-M3, A2-M3 and A3-M2 isomers, which are thought to present the same chemical formula and same molecular weight, have 111.80, 65.20 and 25.88 kJ mol⁻¹ less internal energy than the original asphaltene isomers. These results encourage the development of models based on more stable structures for molecular simulations and emphasize the importance of considering the distribution of π-sextet during this development process.

Table 1 Relative energies^a of the isomers in kJ mol⁻¹ calculated at B3LYP/6-31+G* level of theory

Model	Isomer	ΔE
a The energy is relative to the most stable isomer.
Reduced	A1-M1	111.80
	A1-M2	12.22
	A1-M3	0.00
	A2-M1	65.20
	A2-M2	3.18
	A2-M3	0.00
	A3-M1	25.88
	A3-M2	0.00
Complete	A3-M1	26.51
	A3-M2-a	64.47
	A3-M2-b	0.00

---

Finally, we study the complete asphaltene molecules including the aliphatic chains that were eliminated to simplify the model in the first step of this process. Doing so, A1-M3 and A2-M3 do not suggest any further improvement in the chemical structure, but A3-M2 asphaltene undergoes pentane effect when the aliphatic chains are included. This pentane effect was absent in the simplified model but it appears in the complete one due to the proximity of two aliphatic chains after the new arrangement of aromatic rings in the FAR region. It imposes a need for reconsidering the location of these chains in A3-M2, but keeping the new arrangement of the polycyclic system designed in the previous step. Fig. 6 shows the new isomer in which one aliphatic chain has been moved from position 1 to position 2. Fig. 7, shows A3-M3-a and A3-M3-b new isomers, where the curvature imposed by the pentane effect in A3-M2-a disappears when one of the chains is relocated to position 2. The Energy change due to the relocation of the aliphatic chain is shown in the bottom part of Table 1, where the relative energy of the three complete isomers (including original A3-M1) for the A3 asphaltene is provided. In this case, although the number of π-sextets is increased going from A3-M1 to A3-M2 isomer, the presence of the pentane effect leads to higher energy. This issue is addressed by relocating one of the aliphatic chains leading to 26 kJ mol⁻¹ reduction in internal energy compared to that of the original molecule. Therefore after redesigning the polycyclic aromatic core of asphaltenes using simplified models, it is necessary to check that new steric effect does not appear when the aliphatic chains are included. If the molecules undergo pentane effect with the new arrangement of rings, then the aliphatic chains should be relocated.

Fig. 6 Redesign of A3-M2. (a) A3-M2-a isomer with the original aliphatic chain distribution of A3-M1. (b) New A3-M2-b isomer where one aliphatic chain has been moved from position 1 to position 2, avoiding pentane effect.

Fig. 7 Front and side view of the optimized structures of complete A3 asphaltene isomers at B3LYP/6-31+G* level of theory. The energy relative to the most stable one (A3-M2-b) is shown in parenthesis (in kJ mol⁻¹).

Based on aforementioned results, the new isomers proposed in this paper are shown in Fig. 8. The new geometrical features can be found by comparing their chemical structure to the original asphaltene molecules shown in Fig. 2.

Fig. 8 Asphaltene new molecular structures proposed in the present study. (a) New asphaltene molecules based on previous Greenfield models. (b) and (c) New asphaltene molecules based on Greenfield-modified asphaltenes from Yen–Mullins model.

4. Conclusions

Bringing together DFT calculations and basic organic chemistry concepts of Clar sextet theory, in this paper we propose new asphaltene molecules for molecular simulations of asphalt. Doing so, we provide very simple but elegant rules for a rational design of the arrangement of aromatic rings in the polycyclic aromatic core of asphaltene molecules. This optimizes the distribution of π electrons in the fused aromatic ring region while minimizing geometrical strain within the molecule. The design strategy is divided in three steps. First, simplified molecules with no aliphatic chains are used for the design of the polycyclic aromatic system, which are polycyclic aromatic hydrocarbons that present the same arrangement of aromatic rings as those in asphaltenes. This strategy reduces computational cost for DFT calculations and simplifies the chemical structure under study. Second, the arrangement of rings in polycyclic aromatic core is designed following Clar sextet theory requirements in such a way that the number of Clar π-sextets in the aromatic system is maximized and the new isomers do not undergo pentane effect. Third, the aliphatic chains that were removed in the simplified model are placed back in the rings in their original positions. If the location of the aliphatic chains in the new arrangement of aromatic rings does not undergo any pentane effect, the final isomer is the most stable one. This is the case for A1 and A3 asphaltenes in the present study. On the contrary, if as a consequence of the new arrangement of aromatic rings the aliphatic molecule undergoes pentane effect, a further design step is needed to relocate the chains to address the pentane effect, as it was shown for A3 asphaltene.

Using this approach, isomers that are more stable than the asphaltene molecules used in current asphalt models are proposed. These new isomers (A1-M3, A2-M3 and A3-M1-b) are tens of kJ mol⁻¹ lower in energy than previous asphaltene molecules. This reduction of internal energy was achieved by simultaneously maximizing the number of Clar π-sextets and avoiding the pentane effect. It should be noted that the higher stability of the polycyclic aromatic system of these new molecules is crucial to reproduce the π–π interactions that control the formation of clusters at the nanoscale. We suggest including these new lower-energy isomers in the current asphalt models while considering aforementioned relevant guidelines as a general procedure in future molecular modeling of asphaltenes.

A better understanding of the chemical structure of this class of materials is critical since the molecular interactions at the nanoscale strongly affect rheological and mechanical properties of large-volume applications, including asphalt pavement performance, cracking resistance and ultimately, the durability of infrastructure.

References

1. O. C. Mullins, E. Y. Sheu, A. Hammami and A. G. Marshall, Asphaltenes, Heavy Oils, and Petroleomics, Springer, 2007.
2. D. Lesueur, Adv. Colloid Interface Sci., 2009, 145, 42.
3. O. C. Mullins, Energy Fuels, 2010, 24, 2179.
4. D. D. Li and M. L. Greenfield, Fuel, 2014, 115, 347.
5. O. C. Mullins, H. Sabbah, J. Eyssautier, A. E. Pomerantz, L. Barré, A. B. Andrews, Y. Ruiz-Morales, F. Mostowfi, R. McFarlane, L. Goual, R. Lepkowicz, T. Cooper, J. Orbulescu, R. M. Leblanc, J. Edwards and R. N. Zare, Energy Fuels, 2012, 26, 3986.
6. A. A. Herod, K. D. Bartle, T. J. Morgan and R. Kandiyoti, Chem. Rev., 2012, 112, 3892.
7. I. A. Wiehe, Energy Fuels, 1994, 8, 536.
8. I. A. Wiehe, Energy Fuels, 2012, 26, 4004.
9. Y. Ruiz-Morales and O. C. Mullins, Energy Fuels, 2006, 21, 256.
10. M. R. Gray, Energy Fuels, 2003, 17, 1566.
11. P. Redelius, Energy Fuels, 2004, 18, 1087.
12. L. Zhang and M. L. Greenfield, Energy Fuels, 2007, 21, 1102.
13. L. Zhang and M. L. Greenfield, J. Chem. Phys., 2010, 132, 184502.
14. L. Zhang and M. L. Greenfield, J. Chem. Phys., 2007, 127, 194502.
15. H. Groenzin and O. C. Mullins, Energy Fuels, 2000, 14, 677.
16. I. Kowalewski, M. Vandenbroucke, A. Y. Huc, M. J. Taylor and J. L. Faulon, Energy Fuels, 1996, 10, 97.
17. L. Artok, Y. Su, Y. Hirose, M. Hosokawa, S. Murata and M. Nomura, Energy Fuels, 1999, 13, 287.
18. J. P. Dickie and T. F. Yen, Anal. Chem., 1967, 39, 1847.
19. D. D. Li and M. L. Greenfield, Energy Fuels, 2011, 25, 3698.
20. P. J. Flory, Statistical Mechanics of Chain Molecules, Wiley, New York, 1969.
21. E. Clar, The Aromatic Sextet, Wiley, New York, 1972.
22. M. Randić, Chem. Rev., 2003, 103, 3449.
23. M. Solà, Front. Chem., 2013, 1, 22.
24. Y. Ruiz-Morales, J. Phys. Chem. A, 2004, 108, 10873.
25. Y. Ruiz-Morales and O. C. Mullins, Energy Fuels, 2013, 27, 5017.
26. Y. Ruiz-Morales, Can. J. Chem., 2009, 87, 1280.
27. M. J. T. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery Jr, J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, Ö. Farkas, J. B. Foresman, J. V. OrtizA, J. Cioslowski and D. J. Fox, Revision A.1, Gaussian, Inc., Wallingford, CT, 2009.
28. A. D. Becke, J. Chem. Phys., 1993, 98, 5648.
29. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785.
30. D. Cremer and J. A. Pople, J. Am. Chem. Soc., 1975, 97, 1354.
31. S. Cross, M. M. Kuttel, J. E. Stone and J. E. Gain, J. Mol. Graphics Modell., 2009, 28, 131.
32. M. Kuttel, J. Gain, A. Burger and I. Eborn, J. Mol. Graphics Modell., 2006, 25, 380.
33. W. Humphrey, A. Dalke and K. Schulten, J. Mol. Graphics, 1996, 14, 33.
34. L. Zhang and M. L. Greenfield, Energy Fuels, 2007, 21, 1712.
35. A. Bhasin, R. Bommavaram, M. Greenfield and D. Little, J. Mater. Civ. Eng., 2011, 23, 485.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/c4ra05694a

---