LiquidâVapor Phase Behavior of Asphaltene-like Molecules

Subscriber access provided by University of Oklahoma

Article

Liquid-vapor phase behavior of asphaltene-like molecules Kaustubh Sunil Rane, Vaibhaw Kumar, Scott Wierzchowski, Majeed S. Shaik, and Jeffrey R Errington Ind. Eng. Chem. Res., Just Accepted Manuscript • DOI: 10.1021/ie5035144 • Publication Date (Web): 24 Oct 2014 Downloaded from http://pubs.acs.org on October 28, 2014

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Industrial & Engineering Chemistry Research is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Industrial & Engineering Chemistry Research

Liquid-vapor phase behavior of asphaltene-like molecules by Kaustubh S. Rane1, Vaibhaw Kumar1, Scott Wierzchowski2, * Majeed Shaik3, and Jeffrey R. Errington1, 1

Department of Chemical and Biological Engineering, University at Buffalo, The State University of New York, Buffalo, New York 14260, USA 2

Shell International Exploration and Production Inc., Houston, Texas 77079, USA 3

Shell India Markets Private Limited, Bangalore 560048, India

Submit to: Ind. Eng. Chem. Res.

*

Corresponding author, electronic address: [email protected]

ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules

Abstract Molecular simulation is used to predict liquid-vapor saturation properties of asphaltene-like molecules. We examine nine compounds that are adapted from a structure produced by Boek and coworkers [Energy & Fuels 23, 1209 (2009)] via an optimization technique that generates quantitative molecular representations of asphaltenes. We first consider a sequence of four fused ring structures differentiated by their aromaticity. Next, a collection of structures with a single aromatic core decorated with side chains of variable length and position is examined. Finally, we consider a group of structures with heteroatoms. The General AMBER Force Field (GAFF) and the CHARMM General Force Field (CGenFF) are employed to describe molecular interactions. A general histogram-based approach introduced by Rane and coworkers [J. Chem. Theory Comput. 9, 2552 (2013)] is used to calculate saturated densities, vapor pressures, and enthalpies of vaporization over a wide range of temperatures. Simulation results are compared to estimates from select group contribution methods. Our results indicate that saturated liquid densities and critical temperatures of fused ring structures generally decrease with decreasing aromaticity. The length and distribution of side chains attached to an aromatic core have a significant impact on the saturated liquid densities and a minor impact on the boiling point temperature.

2 ACS Paragon Plus Environment

Page 2 of 45

Page 3 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



I. Introduction Asphaltenes have a significant effect on the thermodynamic and transport properties of crude oil.1,2 These molecules are known to facilitate the clogging of pipes as well as the pores within rocks in oil fields, thereby impeding oil extraction. The identification and design of solutions to mitigate the impact of asphaltenes would benefit from improved modeling capabilities for describing their thermophysical properties. Equations of state currently represent the primary modeling tool employed within the oil and gas industry. This approach relies upon accurate estimates of select liquid-vapor saturation properties for molecules of interest, including the critical temperature, critical pressure, and normal boiling point. Several factors limit the extent to which experiment can be used to evaluate these properties for asphaltenes. Isolating the molecules from solution and identifying their structure has proven difficult. Moreover, the compounds are typically solids at room temperature and thermally breakdown at temperatures well below the boiling and critical temperatures. In the absence of experimental data, alternative approaches are required to estimate key saturation properties. Engineering correlations and group contribution methods are often employed to provide these estimates.3 These approaches are relatively straightforward to implement and enable one to arrive at an estimate quickly with little computational cost. The primary concern when applying this approach to asphaltenes stems from the fact that the molecules of interest have characteristics that typically lie outside of the range of the data set used to fit the parameters. Molecular simulation provides an alternative means to predict saturation properties. Within this approach, one obtains properties by simulating many-body systems governed by a detailed description of the intra- and inter-molecular interactions.4 Effects such as cooperativity and molecular-level organization are captured. As a result, the approach is more fundamental in 3 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules nature, and one expects the range of validity to exceed that of correlation-based methods. However, such calculations involve significant computational resources and often require weeks to complete. Moreover, completing phase equilibria calculations with complex molecules such as asphaltenes is a challenging task. In addition to requiring significant computational resources for computing interactions involving atomistically-detailed versions of these molecules, one must ensure that the simulation algorithms sufficiently sample the intra- and inter-molecular degrees of freedom. Indeed, such simulations push the limits of current methodologies. The molecular structure of asphaltenes has been debated for some time. However, as Mullins details in a recent review,1 consensus appears to be building regarding their size and chemistry. He suggests asphaltenes are characterized by a molecular weight in the 400 to 1000 Da range and are dominated by polycyclic aromatic hydrocarbon (PAH) moieties. Island architectures, which feature a single PAH with 4 to 10 fused rings, are thought to be the most prevalent structures.5 These island cores are typically decorated with short alkyl side chains. Nitrogen, sulfur, and oxygen are also found in small percentages within asphaltenes. Several representative asphaltene molecules and algorithms to produce their quantitative molecular representations have been proposed.1,6-9 Boek et al. recently developed an algorithm that takes as input experimental data, including the molecular weight, elemental analysis, and nuclear magnetic resonance spectroscopy (both 1H and13C data).6 An optimization procedure is used to identify molecular structures that best satisfy these constraints. Boek and coworkers generated several representative structures using this approach, which are consistent with the expectations outlined by Mullins.1 In this work, we use molecular simulation to predict the liquid-vapor saturation properties of several aphaltene-like compounds. We take one of the structures proposed by Boek et al. as a


Page 4 of 45

Page 5 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules basis, and examine how various modifications to this molecule impact the phase behavior. More specifically, we select from the work of Boek et al. a structure with a core PAH consisting of seven fused rings and several side chains. We first focus on the core PAH structure (no side chains) and examine how saturation properties evolve as the aromatic character of the structure is modified. We then investigate how the phase diagram shifts upon modification of the length and distribution of alkyl side chains attached to the PAH core. Finally, we examine the extent to which a nitrogen substitution within the PAH and sulfur substitution within an alkyl side chain impact the phase behavior. One of the key inputs for a classical molecular simulation is the force field, which provides a quantitative description of the intra- and inter-molecular interactions for a collection of molecules. Two commonly used force fields include the General AMBER Force Field (GAFF)10 and the CHARMM General Force Field (CGenFF).11 The strengths and limitations of these models have been explored in several studies.12,13 We recently examined the performance of these models within the context of describing the liquid-vapor saturation properties of naphthalene derivatives (e.g., naphthalene, tetralin, trans-decalin, quinolone, 1-methylnaphthalene), and found that these force fields provide a reasonable description of saturation properties.14 For example, the models captured the critical and boiling temperatures to within average errors of 6 and 4%, respectively. Our results did not indicate that one force field clearly outperforms the other. Here, we use GAFF to perform calculations for all structures and CGenFF for a subset of the compounds. We use molecular simulation to compute saturated densities, vapor pressures, enthalpies of vaporization, and critical properties. We utilize a histogram-based approach15 to obtain coexistence properties over a relatively wide range of temperature. We first complete a series of



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules grand canonical simulations to directly obtain estimates of saturation properties at relatively high temperature. We then use an expanded ensemble procedure that allows us to trace the saturation curve over the entire temperature range of interest. These results are compared with estimates from group contribution methods developed by Rarey and coworkers,16-19 Marrero and Gani,20 and a modified version of the Stein and Brown correlation.21 These group-based methods require as input the structure of the molecule only. These approaches are used to generate boiling temperatures, critical temperatures, pressures, and densities, and the temperature dependence of vapor pressures and enthalpies of vaporization. The paper is organized as follows. In the following section we describe the molecular simulation methods used to compute saturation properties of interest and provide simulation details. We then present the structures we work with and describe the process used to obtain GAFF and CGenFF parameters. Next, we present our simulation results for the molecular models and discuss how various structural factors impact phase behavior. We also discuss the extent to which simulation results and the group contribution methods are consistent. Finally, we provide a few concluding remarks.

II. Molecular Simulation Methods General strategy. We require a means to compute the properties of molecules that possess topologies characterized by a relatively rigid PAH structure, potentially decorated with several flexible alkane side chains (see Figures 1 and 2). The approach that we use has been detailed in earlier reports from our group,14,15 and therefore we limit our discussion here to an overview of the techniques that we employ. We construct phase diagrams via a general two-step histogrambased approach.15 In the first step, we use grand canonical (GC) simulation to locate saturation points at a minimum of two relatively high temperatures. More specifically, GC simulation is 6 ACS Paragon Plus Environment

Page 6 of 45

Page 7 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules used to determine the density probability distribution of a system at a given temperature and activity, and subsequently histogram reweighting22 is employed to locate the coexistence point. This approach provides a so-called “direct” method for determining saturation properties. The aforementioned density probability distribution covers a range of states that span from the homogeneous vapor to homogeneous liquid, including a myriad of inhomogeneous states at intermediate densities.23 This probability distribution typically spans several orders of magnitude and can be difficult to capture. In earlier work,24 we demonstrated that transition matrix Monte Carlo provides an efficient means for computing density probability distributions for “simple” molecules. For molecules with complex topologies and/or interactions, additional challenges come into play. For example, large molecules generally prove difficult to insert and remove from a dense liquid. To address these challenges, we incorporate additional advanced sampling techniques within the histogram-based scheme to arrive at a general method for computing bulk saturation properties.15 Specific strategies that we employ include reservoir grand canonical MC,25 a growth expanded ensemble scheme,26 and hybrid MC methods27 to sample the intramolecular degrees of freedom of molecules. We are interested in collecting coexistence data over a wide range of temperatures. Applying the direct approach becomes computationally demanding and proves inefficient. Moreover, at sufficiently low temperatures the system forms a series of inhomogeneous structures at intermediate densities that prove difficult to sample.23,28 Therefore, we employ an indirect scheme to trace saturation lines over a wide range of temperatures.15 To initiate the calculation, we use coexistence data obtained with the direct method described above to formulate a guess of how properties vary along the liquid-vapor saturation line. More specifically, we now work within an isothermal-isobaric ensemble (NPT), and therefore generate an estimate for how tem-



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules perature and pressure vary along the saturation line by assuming that the logarithm of the vapor pressure ln varies linearly with reciprocal temperature 1/ . States along this path are sampled via isothermal-isobaric expanded ensemble (NPT-TEE) simulations.15 Within these simulations, one creates a series of subensembles that map to specific state points (e.g. 500 subensembles) along the estimate for the saturation line. In this work, we use a scheme in which the subensembles are separated by uniform increments in 1/ . Separate NPT-TEE simulations are completed for the liquid and vapor phases. Thus, we avoid sampling the inhomogeneous states encountered in direct GC simulation at intermediate densities. These simulations provide the relative Gibbs free energy along the path sampled. These data, along with a direct coexistence point, are used to ascertain the accuracy of the linear ln − 1/ guess used to initiate the calculation. As one might suspect, this guess typically proves reasonable at relatively high temperature, but can miss considerably at low temperature. We now modify the relationship. The NPT-TEE simulations provide the density probability distributions sampled in each subensemble. These distributions are used within a histogram reweighting scheme to estimate the change in the Gibbs free energy upon modification of the pressure. The new estimate for at a given temperature is deduced by determining the pressure for which the difference in Gibbs free energy between the liquid and vapor phases is zero. This modified relationship is now used within a new set of NPT-TEE simulations and these data are used to produce estimates of the saturation properties. We typically find that the modified relationship is sufficiently close to the final vapor pressure curve that further simulations are not required. In rare cases, the histogram reweighting process used to estimate the new at a given temperature is pushed outside of its limits and another round of NPT-TEE simulations is required.


Page 8 of 45

Page 9 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules GC simulation details. We employ a cubic simulation cell with an edge length between 48 and 50 Å and span a particle number range from zero to between 100 and 150 molecules. The MC move mix includes molecular displacements and rotations, hybrid MC (HMC) moves (both HMC and HMC-One15), as well as molecular insertions (or growth) and deletions (or reductions). The reservoir includes 500 molecules in an ideal gas state. In previous work, we used a single molecular configuration to populate the reservoir and subsequently equilibrated the molecular configurations before initiating the primary simulation. This scheme typically proves efficient for relatively simple molecules. For the complex molecules studied here, we find that relatively long equilibration times are required. As a result, we adopt a strategy in which we store reservoir configurations at the end of a run, and then use this set of configurations at the outset of the next run. During GC simulations, we regularly complete reservoir MC moves to refresh the reservoir configurations. The total computational cost associated with initializing and maintaining the reservoir typically does not exceed 15% of the total computational effort. The Lo and Palmer potential29 is employed to describe the manner in which a partial molecule interacts with the rest of the system during the growth expanded ensemble process. The number of growth stages employed spans linearly from one at low density to a molecule-specific maximum for dense liquid conditions. For the molecules studied here, the upper bound is set between 20 and 30. A spherical cutoff coupled with a mean-field correction is used to account for long-range Lennard-Jones interactions.4 The Ewald summation is used to compute electrostatic interactions.4 Lennard-Jones interactions and the real space contribution to the Ewald summation are cutoff at a distance of 16 Å. The algorithm described by Frenkel and Smit4 with a precision of 5 × 10 is used to determine the damping factor and number of lattice vectors associated with the Ewald summation.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 45

Liquid-vapor phase behavior of asphaltene-like molecules NPT-TEE simulations. We employ a cubic simulation cell consisting of 100 to 150 molecules. The MC move mix includes molecular displacements and rotations, HMC, and HMC-One moves, volume change attempts, and subensemble change attempts. Volume changes are selected with a probability of approximately 1/2 and subensemble changes are attempted with a probability of 1%. The spherical cutoff distance for the Lennard-Jones and real-space Ewald contribution is taken to be a constant fraction of the simulation cell length, usually 25%. The damping factor and number of lattice vectors are determined using the approach described in the previous paragraph. Simulations are performed over a broad range of = 1/ ( is Boltzmann’s constant) with a discretization of Δ = 1.0 × 10 J–1. We take as the order parameter because the Gibbs free energy scales approximately linearly with inverse temperature along the liquid-vapor saturation curve. Note that the temperatures of successive bins and ! are related as 1/ ! − 1/ = Δ. A typical temperature difference for adjacent subensembles is 1 K

III. Molecular Models Structures. The compounds examined in this study are presented in Figure 1. As noted above, the structures we consider are based upon an asphaltene molecule generated by Boek et al.6 The structure from their work is provided in the lower left corner of Figure 1. We first make a few changes to arrive at the C37H33NS structure shown in the lower right corner of Figure 1. Specifically, we (1) remove the methyl group pointing toward an aromatic hydrogen (right side of structure), (2) shift the lower left methyl group, (3) make the lower left ring aromatic, (4) shorten the long sulfur-containing side chain, and (5) swap the lower right propyl group for an ethyl group. The latter two modifications are implemented for computational efficiency. The first two modifications are motivated by the work of Li and Greenfield.30 They identified several geometric factors that lead to high energy structures. Several of the issues they highlight are related to the 10 ACS Paragon Plus Environment

Page 11 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules “pentane effect”, which is associated with steric repulsion between atoms separated by six bond lengths. The canonical example relates to the interaction between two hydrogens stemming from the terminal carbons of a pentane molecule arranged in a ring-like geometry. For the Boek et al. structure, this effect comes into play when the methyl group on the right side of the molecule conflicts with an aromatic hydrogen. Li and Greenfield also identified geometries with side chains positioned in a manner that is consistent with the two left-most side chains as being relatively high in energy. In their work, they considered the specific case in which two ethyl groups stem from equivalent positions. Again, one can rationalize these high energy configurations by considering the pentane effect. To mitigate this issue, we shift the methyl group by one position on the aromatic ring.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 45


Figure 1. Molecules examined in this study. MW denotes the molecular weight of the molecule.

After arriving at the C37H33NS structure, we modify this molecule in various ways, such that we can systematically study the effect of core aromaticity, the length and positioning of side chains, and the presence of heteroatoms on liquid-vapor saturation properties. We first eliminate the side chains and heteroatoms and subsequently modify the aromatic character of the resulting core. The four structures examined are provided in the first row of Figure 1. We then retain the completely aromatic core and add several side chains (middle row of Figure 1). The C34H26 and C39H36 structures are differentiated by the length of the side chains, whereas the two C39H36 structures are differentiated by the distribution of the side chains. Finally, we consider the presence of a nitrogen substitution (C38H35N) and an additional sulfur substitution (C37H33NS). 12 ACS Paragon Plus Environment

Page 13 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules Quantum mechanical calculations are completed to obtain optimized geometries for each of the structures. Specifically, we examine the molecules in vacuum using B3PW91 level of KohnSham density functional theory (DFT) and the 6-31+G** basis set with the commercially available Q-Chem software (version 4.0).31 Snapshots of C30H18 and C37H33NS are provided in Figure 2. The aromatic core exhibits considerable curvature. The rationale for this curvature is again related to the pentane effect. In three locations within the fused ring structure aromatic hydrogens separated by six bond lengths occupy a similar region of space within the planar configuration. As a result, the rings twist out of plane to relieve the steric repulsion.

Figure 2. Optimized molecular structures for C30H18 and C30H18NS. Gray, white, blue, and yellow atoms represent carbon, hydrogen, nitrogen, and sulfur, respectively.

Force field. We are interested in a force field that provides interaction parameters for fused aromatic ring systems, including those with heteroatoms and alkyl side chains. Moreover, it is de-



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 45

Liquid-vapor phase behavior of asphaltene-like molecules sirable to have a force field that is parameterized to reproduce liquid-vapor saturation properties of molecules with similar chemistries over a wide range of temperature. However, to our knowledge, there is no single force field that satisfies both of these criteria. We considered using the Transferable Potentials for Phase Equilibria Force Field (TraPPE), which is parameterized to accurately capture liquid-vapor saturation properties. While TraPPE has been extended recently to describe the properties of complex molecules like aromatic compounds containing heteroatoms32 and PAHs belonging to Kata-annellated systems (e.g., naphthalene and anthracene),33 it is not clear to us how to extend TraPPE to model the molecules studied here. GAFF was developed as an extension of the Amber force field for drug-like molecules.10 It provides interaction parameters for a wide variety of heterocyclic ring systems and functional groups. Also, the developers have provided a straightforward procedure to obtain parameters related to the force field.34 Bonded parameters and partial charges on atoms are fit to data from quantum calculations and experimental crystal structure data. Dispersion interaction parameters are fit to reproduce liquid densities and enthalpies of vaporization obtained from experiments. The force field has proven to be helpful in studies involving small molecules.14 CGenFF is an extension of the CHARMM force field and, like GAFF, was developed for performing molecular simulations with drug-like molecules. Bonded interaction parameters are fit to the results from quantum calculations. CGenFF dispersion parameters are fit to reproduce liquid densities and enthalpies of vaporization obtained from experiments.11 Like GAFF, CGenFF has been widely used, mainly to study the behavior of small molecules in biological systems. With respect to the current study, the key limitations of GAFF and CGenFF are that they have not been parameterized to liquid-vapor saturation properties over a wide range of temperature. However, this concern is mitigated to some extent by the results from a recent study by our group, wherein we con-


Page 15 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules sider the ability of GAFF and CGenFF to capture the saturation properties of various naphthalene derivatives.14 Our results indicate that GAFF and CGenFF are just slightly inferior to TraPPE in their ability to capture saturation properties. In this study, we use GAFF to predict phase diagrams for all of the structures provided in Figure 1. CGenFF is used to construct phase diagrams for C30H18 and C30H48. GAFF parameters are obtained with the help of the Antechamber program35 and CGenFF interaction parameters are obtained via the Paramchem website.34 The partial atomic charges for GAFF are calculated using AM1-BCC method,36 as provided in Antechamber.

IV. Results and Discussion We begin by illustrating the methods used to construct phase diagrams. First, GC simulation is used to obtain the particle number (density) probability distribution Π at a given temperature and activity. At conditions close to coexistence, this distribution is bimodal, with the two peaks associated with the vapor and liquid phases. The distribution is then shifted using histogram reweighting to an activity for which the areas under the two peaks are equal, which one can show is equivalent to the pressures of the two phases being equal.24 The resulting distributions for C30H26 at temperatures of 900 and 950 K are provided in Figure 3. Thermodynamic properties are extracted from the probability distribution.24 For example, saturated densities are related to the first moments of the peaks and the vapor pressure is related to the probability of observing an empty simulation box (Π = 0).24



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 45


Figure 3. Particle number probability distributions from grand canonical Monte Carlo simulation for C30H26 at saturation conditions.

Results from the GC simulations are used to formulate a guess for the pressure-temperature relationship along the saturation line. It is assumed that the logarithm of the vapor pressure ln is a linear function of the reciprocal temperature 1/ . The GC points and the guess function for C30H26 are provided in Figure 4. Isothermal-isobaric expanded ensemble simulations that span this line are performed for the vapor and liquid phases. These simulations provide the relative Gibbs free energy of the respective phase along this path. One of the two GC points serves as a reference temperature at which the Gibbs free energies of the two phases must be equal. This construction allows one to compute the Gibbs free energy difference at each point along the path. The shift in the pressure required to eliminate the free energy difference is now estimated. A number of techniques can be used to complete this calculation. Formally, one uses volume probability distributions collected within each subensemble to describe the change in the 16 ACS Paragon Plus Environment

Page 17 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules Gibbs free energy with pressure.15 Alternatively, simple thermodynamic relations can be used to estimate this relationship. For example, if one takes the liquid to be incompressible, the free energy difference is given by Δ# = 〈%〉Δ, where Δ# is the extensive change in free energy, is the number of molecules in the expanded ensemble simulation, 〈%〉 is the ensemble average of the molar volume, and Δ is the shift in the pressure. For the vapor, one might take the gas to be ideal, and arrive at a simple relation for the variation in Δ# with . Regardless of the approach used, one now arrives at a new estimate for . A new set of expanded ensemble simulations is completed, providing a means to further refine the estimate. The converged relation for C30H26 is provided in Figure 4.

Figure 4. Schematic of the process used to trace saturation curves. Open circles represent the results from direct grand canonical simulations. The dashed and solid curves represent the guess and converged, respectively, relationships associated with the isothermal-isobaric temperature expanded ensemble approach.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 45

Liquid-vapor phase behavior of asphaltene-like molecules We now turn our attention to the phase behavior of the asphaltene-like molecules. In what follows, we provide simulation results obtained with GAFF and CGenFF. We also present estimates of select saturation properties generated using the group contribution methods of Rarey and coworkers,16-19 Marrero and Gani,20 and Stein and Brown.21 The Rarey and coworkers (NMRR) method includes contributions from individual chemical groups, second-order correction terms (often related to molecular geometry and/or connectivity), and pairwise group interaction terms. The authors provide working relations for estimating the boiling temperature ' ,19 critical temperature ( , pressure ( , and density )( ,18 and temperature dependence of the vapor pressure16,17 (from which the enthalpy of vaporization is obtained). The group has introduced two correlations for the vapor pressure. In this work, we show results for the Moller et al. version only.17 These methods are available in the ARTIST software provided by Dortmund Data Bank Software & Separation Technology GmbH.37 Here, we use ARTIST to identify the groups relevant to a molecule and subsequently compute the properties manually. The Marrero and Gani method20 (MG) includes primary-level contributions from simple groups as well as higherlevel contributions focused on polyfunctional and structural groups related to molecular fragments. They provide working relations for calculating ' , ( , ( , and )( . The parameters for their approach are obtained via consideration of large complex polycyclic molecules, rendering the method particularly relevant to this study. The third group contribution technique considered is the “adapted Stein and Brown method” (aSB) available within the Estimation Program Interface (EPI) Suite (version 4.11) provided by the United States Environmental Protection Agency.38 This approach is based upon the original Stein and Brown method,21 which, in turn, is an extension of the Joback method.39 The adapted version includes additional functional groups and correction factors. In this work, we present estimates for ' , which are obtained via the EPI


Page 19 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules software suite.38 Finally, we show estimates for saturated liquid densities based upon the Reid et al. correlation.40 This equation requires as input ( , )( , and the acentric factor *, which we obtain from the methods of Rarey and coworkers. We first explore how the aromaticity of the fused-ring core structure impacts the phase diagram. Saturated densities ), vapor pressures , and enthalpies of vaporization Δ+,- for C30H18, C30H26, C30H34, and C30H48 are provided in Figure 5. Results for GAFF and CGenFF (C30H18 and C30H48 only) as well as the NMRR (critical properties, , and Δ+,- ) and Reid et al. (). / ) correlations are shown. Values for ' , ( , ( , and )( are also collected in Table 1. As the aromaticity increases, the saturated liquid densities increase in a regular manner. The GAFF, CGenFF, and correlation-based ). / results for the fully aromatic C30H18 are reasonably consistent, with the correlation providing values that deviate by less than 4% from the simulation data at all temperatures. Similar agreement is obtained for C30H26 and C30H34. In contrast, GAFF, CGenFF, and correlation-based results for the fully saturated C30H48 differ considerably. The NMRR estimate for )( decreases appreciably in going from the partially aromatic to fully saturated fused ring system, thereby pulling the ). / curve to much lower values. Within the Reid et al. correlation, ). / scales linearly with )( , so the 30% decrease in )( from C30H34 to C30H48 leads to a similar shift in the ). / curve.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 45


Figure 5. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C30H18, C30H26, C30H34, and C30H48. Solid and dashed curves represent simulation data obtained with GAFF and CGenFF, respectively. Diamond, circle, and × symbols correspond to critical points for GAFF, CGenFF, and the NMRR correlation, respectively. Up and down triangles correspond to vapor pressures and enthalpies of vaporization, respectively, generated with the NMRR correlation. Squares represent predictions from the Reid et al. correlation.


Page 21 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules Table 1. Critical and boiling point properties of the molecules studied in this work. Compound C30H18

C30H26

C30H34

C30H48

C34H26

C39H36 (L)

C39H36 (S)

C38H35N

C37H33NS

Model GAFF CGenFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF CGenFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB

21 (kg/m3) 3342 3254 367 299 --3104 349 295 --2861 333 291 --2662 2563 293 220 --3012 346 274 --2792 328 253 --2766 328 253 --2753 321 252 --3515 333 259 ---

01 (K) 10984 11339 1125 1122 --10719 1085 1105 --10474 1038 1088 --11058 9449 918 914 --10391 1135 1104 --10298 1146 1098 --10057 1145 1085 --106611 1143 1105 --10237 1149 1134 ---

31 (MPa) 2.51 2.72 1.5 2.1 --2.12 1.5 1.8 --2.10 1.5 1.6 --1.31 2.14 1.3 1.1 --1.60 1.2 1.5 --1.41 1.1 1.1 --1.21 1.1 1.0 --1.63 1.2 1.1 --0.80 1.2 1.3 ---

04 (K) 8390 8532 833 866 885 8221 801 852 849 8152 763 837 813 8202 7897 704 714 730 8341 848 893 931 8463 865 924 989 8321 866 921 989 8712 857 930 995 8821 862 947 1017

GAFF results suggest that ( and ' exhibit a minimum upon varying the aromaticity. While CGenFF simulations are completed for the extreme cases only, these results suggest that ( and ' decrease with decreasing aromaticity. All of the group contribution methods suggest trends that are qualitatively consistent with CGenFF. The different trends captured by GAFF and CGenFF regarding the influence of aromaticity on ( and ' are consistent with our work 21 ACS Paragon Plus Environment


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 45

Liquid-vapor phase behavior of asphaltene-like molecules focused on naphthalene derivatives.14 In this earlier study, we examined naphthalene, tetralin, and trans-decalin, molecules consisting of two, one, and zero aromatic rings, respectively. For these compounds, experimental data and CGenFF indicate that ( and ' decrease with decreasing aromaticity. In contrast, GAFF suggests that these temperatures exhibit a minimum upon variation of the aromaticity. Collectively, these results suggest that GAFF does not accurately differentiate between aromatic and saturated fused ring compounds. GAFF vapor pressure curves for the C30H18, C30H26, and C30H34 fused ring structures are quantitatively similar over a fairly wide temperature range. For temperatures above 700 K, values do not differ by more than a factor of three. As a result, the enthalpies of vaporization are also fairly consistent, particularly at moderate to low temperature. The GAFF and CGenFF descriptions of C30H18 are relatively consistent, with nearly identical and Δ+,- curves. In contrast, GAFF and CGenFF and Δ+,- curves for C30H48 differ significantly. The curves are offset at relatively high temperature and the CGenFF data suggest a sharper decrease in with decreasing temperature. The latter trend results in CGenFF Δ+,- values that are approximately 30% larger than those for GAFF. The NMRR correlation provides estimates for the partially aromatic ring structures that are reasonably consistent with GAFF data. Similar to GAFF, the correlation suggests that there is relatively little variation in and Δ+,- curves for C30H18, C30H26, and C30H34. Correlation data for C30H48 are noticeably displaced from those for the other compounds. We briefly comment on the structure of the fluids examined here. Asphaltenes are known to form nanoaggregates, with the PAH cores of the molecules loosely stacking together.1 Scattering experiments have provided evidence for such structures in solution.41,42 For the conditions explored in this work, one might expect to observe nanoaggregates in the vapor phase. Also, 22 ACS Paragon Plus Environment

Page 23 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules given the conventional view of PAHs as planar structures, one might expect the molecules within the liquid to organize into a nematic liquid crystal at sufficiently low temperature. However, we do not observe either of these assemblies here. The absence of such structuring is likely related to the twisted nature of the PAH core, which restricts the ability of the molecules to assemble into organized structures. We now examine the influence of alkyl side chains on the phase behavior of asphaltene structures. Figure 6 contains saturation properties for C30H18, C34H26, C39H36 (S), and C39H36 (L). These molecules are differentiated by the length and/or distribution of alkyl chains attached to a fully aromatic PAH core. As the total carbon number associated with the alkyl chains increases, the saturated liquid densities decrease. This result is reasonable, as one finds the densities approaching those characteristic of saturated hydrocarbons with increasing alkyl chain length. The one case examined here to study the role of chain distribution (C39H36 (S) versus C39H36 (L)) suggests that varying the chain distribution at constant carbon number has relatively little impact on the saturated liquid density profile. The correlation-based description of the impact of alkyl side chains on ). / is qualitatively consistent with GAFF. Interestingly, the critical temperature decreases as the total chain carbon number increases. Presumably, the alkyl chains frustrate the ability of the PAH cores to associate, thereby reducing the overall strength of the intermolecular interaction. The NMRR correlation suggests a trend that is qualitatively consistent with GAFF, while the MG correlation suggests the opposite (see Table 1). As noted above, the aromatic PAH core examined here has a pronounced twisted geometry. It would be interesting to examine the role of alkyl chains with a relatively flat PAH core structure that promotes core stacking.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 45


Figure 6. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C30H18, C34H26, C39H36 (S), and C39H36 (L). Curves and symbols are defined in the same manner as in Figure 5.

Similar to the group of molecules with varying aromaticity examined above, we find that the GAFF vapor pressure curves are quantitatively similar. The boiling temperatures fall within the relatively narrow range of 832 to 846 K and the curves exhibit a similar shape. The enthalpy of vaporization curves are also relatively consistent, with no clear trend related to alkyl chain 24 ACS Paragon Plus Environment

Page 25 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules length or distribution. We note that in our previous work focused on naphthalene derivatives, we examined 1-methylnaphthalene, which consists of a methyl group attached to the naphthalene ring structure. In this case, GAFF accurately captures the shift in saturation properties related to the addition of the methyl group, lending some optimism that the model is providing reasonable results for the decorated PAH structures examined here. Relative to GAFF, the NMRR correlation suggests that alkyl chains have a more pronounced impact on ' . More specifically, the correlation predicts a 60 K increase in ' as a result of adding alkyl chains to the aromatic core. However, GAFF and the NMRR correlation predict similar temperature dependencies for , resulting in a relatively similar description of Δ+,- . Asphaltene molecules often contain nitrogen, sulfur, and oxygen, which appear as substitutes for carbon atoms within the ring structure or side chains. Figure 7 provides saturation properties for two molecules containing nitrogen and/or sulfur substitutions, C38H35N and C37H33NS, along with C39H36 (L), which possesses the same molecular architecture as the others. Our results indicate that these substitutions lead to an increase in the normal boiling temperature and enthalpies of vaporization, suggesting that the substitutions result in a stronger average interaction between molecules. The saturated liquid densities also increase as a result of the nitrogen and sulfur substitutions. The NMRR correlation provides trends that are qualitatively consistent with GAFF simulation results.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 45


Figure 7. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C39H36 (L), C38H35N, and C37H33NS. Curves and symbols are defined in the same manner as in Figure 5.

We now consider further the extent to which molecular simulation and group contribution methods provide a consistent description of the saturation properties of asphaltene-like molecules. We highlight the normal boiling temperature, as this is a focal point for many of the group contribution methods. Figure 8 contains estimates of ' for the compounds studied here. All of 26 ACS Paragon Plus Environment

Page 27 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules the correlations suggest that ' monotonically decreases with the aromaticity of the fused ring structure, with estimates for the total shift in ' spanning from 129 to 155 K. Analogous shifts for GAFF and CGenFF are 19 and 64 K, respectively. The MG correlation suggests that ' scales approximately linearly with the number of aromatic rings in the fused ring structure, while the NMRR and aSB correlations predict a slower initial decline in ' with decreasing aromatic character followed by a sharp decrease when shifting to the fully saturated system. With the exception of C30H48, the GAFF trend is most consistent with the NMRR correlation. There is more diversity in the results related to alkyl side chains. While all of the correlations suggest that the addition of alkyl side chains leads to an increase in ' , the magnitude of the estimate varies considerably, with the MG, NMRR, and aSB methods predicting increases in ' of 32, 58, and 104 K, respectively, in moving from the C30H18 to C39H36 (L) structure. As noted above, GAFF predicts that the alkyl side chains have relatively little influence on ' . This difference could be related to a packing effect that is captured by the molecular simulations but not the group contribution methods. The group contribution methods also show diversity with respect to the influence of the nitrogen and sulfur substitutions. The MG, NMRR, and aSB methods predict shifts in ' of –3, 23, and 28 K, respectively, in moving from the C39H36 (L) to C37H33NS structure. GAFF provides an analogous shift in ' of 36 K. These results indicate that the correlations underestimate relative to GAFF the enhancement in the strength of the intermolecular interaction due to nitrogen and sulfur substitutions. Interestingly, Zhang and Greenfield reported a relatively small contribution to the polar interaction associated with nitrogen- and sulfur-containing molecules in a model asphalt mixture from analysis of Hansen solubility parameter correlations.43



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 45


Figure 8. Normal boiling temperature. For a given molecule, the bars from left to right represent GAFF (solid blue), CGenFF (red diagonal hatch, not present for all molecules), MG (green square hatch), NMRR (violet horizontal hatch), aSB (maroon checkerboard).


Page 29 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



V. Conclusions We have employed molecular simulation to study the liquid-vapor phase behavior of asphaltene-like molecules. We examined nine compounds that are based upon structures generated by Boek and coworkers.6 These compounds were selected such that we could study the impacts of fused ring aromaticity, alkyl side chains, and heteroatoms on saturation properties. Molecular simulation was used to construct phase diagrams for all compounds with interactions described by GAFF and a subset of the compounds was examined with CGenFF. A general histogrambased approach was used to calculate saturated densities, vapor pressures, and enthalpies of vaporization. Simulation results were compared to results from the group contribution methods of Rarey and coworkers,16-19 Marrero and Gani,20 and Stein and Brown.21 The results generated here suggest that the saturated liquid density, boiling temperature, and critical temperature increase with increasing aromaticity. The force fields and group contribution methods examined provided a reasonably consistent description of the fully aromatic fused ring structure, but differed considerably for the fully saturated ring structure. The addition of alkyl side chains to the PAH core leads to a decrease in critical temperature and saturated liquid density. Nitrogen and sulfur substitutions resulted in an increase in the normal boiling temperature. The results generated here suggest that molecular simulation can serve as a powerful tool for understanding the thermophysical properties of asphaltenes. The methods employed here could be used to study the phase behavior of a more diverse set of structures (e.g. other “island” as well as “archipelago” architectures). Modified versions of these methods could be employed to study more complex phase behavior, such as the precipitation of asphaltenes from solution. Molecular simulation could also prove useful in obtaining a better understanding of the impact of



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules asphaltenes on wetting and other interfacial phenomena. For example, one could probe the extent to which surface active asphaltenes modify the manner in which water and oil wet mineral surfaces.

Acknowledgements We gratefully acknowledge financial support of Shell India Markets Private Limited. Computational resources were provided in part by the University at Buffalo Center for Computational Research and the Rensselaer Polytechnic Institute Computational Center for Nanotechnology Innovations.

Associated Content Supporting Information The data presented in Figures 5-8 are provided. This material is available free of charge via the Internet at http://pubs.acs.org.


Page 30 of 45

Page 31 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



References (1) Mullins, O. C. The Asphaltenes. Annu. Rev. Anal. Chem. 2011, 4, 393-418. (2) C&En Cover Story: Digging into Asphaltenes http://pubs.acs.org/cen/coverstory/87/8738cover.html (accessed on 04/22/2014). (3) van Speybroeck, V.; Gani, R.; Meier, R. J. The Calculation of Thermodynamic Properties of Molecules. Chem. Soc. Rev. 2010, 39, 1764-1779. (4) Frenkel, D.; Smit, B. Understanding Molecular Simulation; 2nd ed.; Academic Press: London, 2001; Vol. 1. (5) Ruiz-Morales, Y.; Mullins, O. C. Polycyclic Aromatic Hydrocarbons of Asphaltenes Analyzed by Molecular Orbital Calculations with Optical Spectroscopy. Energy & Fuels 2006, 21, 256-265. (6) Boek, E. S.; Yakovlev, D. S.; Headen, T. F. Quantitative Molecular Representation of Asphaltenes and Molecular Dynamics Simulation of Their Aggregation†. Energy & Fuels 2009, 23, 1209-1219. (7) Neurock, M.; Nigam, A.; Trauth, D.; Klein, M. T. Molecular Representation of Complex Hydrocarbon Feedstocks through Efficient Characterization and Stochastic Algorithms. Chemical Engineering Science 1994, 49, 4153-4177. (8) Kowalewski, I.; Vandenbroucke, M.; Huc, A. Y.; Taylor, M. J.; Faulon, J. L. Preliminary Results on Molecular Modeling of Asphaltenes Using Structure Elucidation Programs in Conjunction with Molecular Simulation Programs. Energy & Fuels 1996, 10, 97107. (9) Sheremata, J. M.; Gray, M. R.; Dettman, H. D.; McCaffrey, W. C. Quantitative Molecular Representation and Sequential Optimization of Athabasca Asphaltenes. Energy & Fuels 2004, 18, 1377-1384. (10) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General Amber Force Field. J. Comput. Chem. 2004, 25, 1157-1174. (11) Vanommeslaeghe, K.; Hatcher, E.; Acharya, C.; Kundu, S.; Zhong, S.; Shim, J.; Darian, E.; Guvench, O.; Lopes, P.; Vorobyov, I.; Mackerell, A. D. Charmm General Force Field: A Force Field for Drug-Like Molecules Compatible with the Charmm All-Atom Additive Biological Force Fields. J. Comput. Chem. 2010, 31, 671-690. (12) Mobley, D. L.; Bayly, C. I.; Cooper, M. D.; Shirts, M. R.; Dill, K. A. Small Molecule Hydration Free Energies in Explicit Solvent: An Extensive Test of Fixed-Charge Atomistic Simulations. J. Chem. Theory Comput. 2009, 5, 350-358. (13) Caleman, C.; van Maaren, P. J.; Hong, M.; Hub, J. S.; Costa, L. T.; van der Spoel, D. Force Field Benchmark of Organic Liquids: Density, Enthalpy of Vaporization, Heat Capacities, Surface Tension, Isothermal Compressibility, Volumetric Expansion Coefficient, and Dielectric Constant. J. Chem. Theory Comput. 2011, 8, 61-74. (14) Kumar, V.; Rane, K. S.; Wierzchowski, S.; Shaik, M.; Errington, J. R. Evaluation of the Performance of GAFF and CGenFF in the Prediction of Liquid–Vapor Saturation Properties of Naphthalene Derivatives. Ind. Eng. Chem. Res. 2014, 53, 16072-16081. (15) Rane, K. S.; Murali, S.; Errington, J. R. Monte Carlo Simulation Methods for Computing Liquid–Vapor Saturation Properties of Model Systems. J. Chem. Theory Comput. 2013, 9, 2552-2566.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Liquid-vapor phase behavior of asphaltene-like molecules (16) Nannoolal, Y.; Rarey, J.; Ramjugernath, D. Estimation of Pure Component Properties: Part 3. Estimation of the Vapor Pressure of Non-Electrolyte Organic Compounds Via Group Contributions and Group Interactions. Fluid Phase Equilibr. 2008, 269, 117-133. (17) Moller, B.; Rarey, J.; Ramjugernath, D. Estimation of the Vapour Pressure of Non-Electrolyte Organic Compounds Via Group Contributions and Group Interactions. J. Mol. Liq. 2008, 143, 52-63. (18) Nannoolal, Y.; Rarey, J.; Ramjugernath, D. Estimation of Pure Component Properties: Part 2. Estimation of Critical Property Data by Group Contribution. Fluid Phase Equilibr. 2007, 252, 1-27. (19) Nannoolal, Y.; Rarey, J.; Ramjugernath, D.; Cordes, W. Estimation of Pure Component Properties: Part 1. Estimation of the Normal Boiling Point of Non-Electrolyte Organic Compounds Via Group Contributions and Group Interactions. Fluid Phase Equilibr. 2004, 226, 45-63. (20) Marrero, J.; Gani, R. Group-Contribution Based Estimation of Pure Component Properties. Fluid Phase Equilibr. 2001, 183–184, 183-208. (21) Stein, S. E.; Brown, R. L. Estimation of Normal Boiling Points from Group Contributions. J. Chem. Inf. Comput. Sci. 1994, 34, 581-587. (22) Ferrenberg, A. M.; Swendsen, R. H. New Monte Carlo Technique for Studying Phase Transitions. Phys. Rev. Lett. 1988, 61, 2635-2638. (23) MacDowell, L. G.; Shen, V. K.; Errington, J. R. Nucleation and Cavitation of Spherical, Cylindrical, and Slablike Droplets and Bubbles in Small Systems. J. Chem. Phys. 2006, 125, 034705-034715. (24) Errington, J. R. Direct Calculation of Liquid–Vapor Phase Equilibria from Transition Matrix Monte Carlo Simulation. J. Chem. Phys. 2003, 118, 9915-9925. (25) Errington, J. R.; Panagiotopoulos, A. Z. New Intermolecular Potential Models for Benzene and Cyclohexane. J. Chem. Phys. 1999, 111, 9731-9738. (26) Escobedo, F. A.; de Pablo, J. J. Expanded Grand Canonical and Gibbs Ensemble Monte Carlo Simulation of Polymers. J. Chem. Phys. 1996, 105, 4391-4394. (27) Duane, S.; Kennedy, A. D.; Pendleton, B. J.; Roweth, D. Hybrid Monte Carlo. Phys. Lett. B 1987, 195, 216-222. (28) Godawat, R.; Jamadagni, S. N.; Errington, J. R.; Garde, S. Structure, Stability, and Rupture of Free and Supported Liquid Films and Assemblies in Molecular Simulations. Ind. Eng. Chem. Res. 2008, 47, 3582-3590. (29) Lo, C.; Palmer, B. Alternative Hamiltonian for Molecular Dynamics Simulations in the Grand Canonical Ensemble. J. Chem. Phys. 1995, 102, 925-931. (30) Li, D. D.; Greenfield, M. L. High Internal Energies of Proposed Asphaltene Structures. Energy Fuels 2011, 25, 3698-3705. (31) Shao, Y.; Molnar, L. F.; Jung, Y.; Kussmann, J.; Ochsenfeld, C.; Brown, S. T.; Gilbert, A. T. B.; Slipchenko, L. V.; Levchenko, S. V.; O'Neill, D. P.; DiStasio Jr, R. A.; Lochan, R. C.; Wang, T.; Beran, G. J. O.; Besley, N. A.; Herbert, J. M.; Yeh Lin, C.; Van Voorhis, T.; Hung Chien, S.; Sodt, A.; Steele, R. P.; Rassolov, V. A.; Maslen, P. E.; Korambath, P. P.; Adamson, R. D.; Austin, B.; Baker, J.; Byrd, E. F. C.; Dachsel, H.; Doerksen, R. J.; Dreuw, A.; Dunietz, B. D.; Dutoi, A. D.; Furlani, T. R.; Gwaltney, S. R.; Heyden, A.; Hirata, S.; Hsu, C.-P.; Kedziora, G.; Khalliulin, R. Z.; Klunzinger, P.; Lee, A. M.; Lee, M. S.; Liang, W.; Lotan, I.; Nair, N.; Peters, B.; Proynov, E. I.; Pieniazek, P. A.; Min Rhee, Y.; Ritchie, J.; Rosta, E.; David Sherrill, C.; Simmonett, A. C.; Subotnik, J. E.; Lee Woodcock Iii, H.; Zhang, W.; Bell, 32 ACS Paragon Plus Environment

Page 32 of 45

Page 33 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Liquid-vapor phase behavior of asphaltene-like molecules A. T.; Chakraborty, A. K.; Chipman, D. M.; Keil, F. J.; Warshel, A.; Hehre, W. J.; Schaefer Iii, H. F.; Kong, J.; Krylov, A. I.; Gill, P. M. W.; Head-Gordon, M. Advances in Methods and Algorithms in a Modern Quantum Chemistry Program Package. Phys. Chem. Chem. Phys. 2006, 8, 3172-3191. (32) Rai, N.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 9. Explicit Hydrogen Description of Benzene and Five-Membered and Six-Membered Heterocyclic Aromatic Compounds. J. Phys. Chem. B 2007, 111, 10790-10799. (33) Rai, N.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 10. ExplicitHydrogen Description of Substituted Benzenes and Polycyclic Aromatic Compounds. J. Phys. Chem. B 2012, 117, 273-288. (34) Paramchem. http://www.paramchem.org (accessed on 01/29/2014). (35) Wang, J.; Wang, W.; Kollmann, P.; Case, D. Antechamber, an Accessory Software Packagefor Molecular Mechanical Calculation. J. Comput. Chem. 2005, 25, 11571174. (36) Jakalian, A.; Bush, B. L.; Jack, D. B.; Bayly, C. I. Fast, Efficient Generation of High-Quality Atomic Charges. Am1-Bcc Model: I. Method. J. Comput. Chem. 2000, 21, 132146. (37) Artist. http://www.ddbst.com/artist-property-estimation.html ( accessed on 04/18/2014). (38) Us Epa. [2014]. Estimation Programs Interface Suitetm http://www.epa.gov/oppt/exposure/pubs/episuite.htm (accessed on 4/18/2014). (39) Joback, K. G.; Reid, R. C. Estimation of Pure-Component Properties from GroupContributions. Chem. Eng. Commun. 1987, 57, 233-243. (40) Reid R. C., P. J. M., Sherwood T. K. The Properties of Gases and Liquids; 5 ed.; McGraw-Hill: New York, 1977. (41) Gawrys, K. L.; Blankenship, G. A.; Kilpatrick, P. K. Solvent Entrainment in and Flocculation of Asphaltenic Aggregates Probed by Small-Angle Neutron Scattering. Langmuir 2006, 22, 4487-4497. (42) Gawrys, K. L.; Kilpatrick, P. K. Asphaltenic Aggregates Are Polydisperse Oblate Cylinders. J. Colloid Interface Sci. 2005, 288, 325-334. (43) Zhang, L.; Greenfield, M. L. Effects of Polymer Modification on Properties and Microstructure of Model Asphalt Systems. Energy & Fuels 2008, 22, 3363-3375.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 45


Tables Table 1. Critical and boiling point properties of the molecules studied in this work. Compound C30H18

C30H26

C30H34

C30H48

C34H26

C39H36 (L)

C39H36 (S)

C38H35N

C37H33NS

Model GAFF CGenFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF CGenFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB GAFF MG NMRR aSB

01 (K) 10984 11339 1125 1122 --10719 1085 1105 --10474 1038 1088 --11058 9449 918 914 --10391 1135 1104 --10298 1146 1098 --10057 1145 1085 --106611 1143 1105 --10237 1149 1134 ---

21 (kg/m3) 3342 3254 367 299 --3104 349 295 --2861 333 291 --2662 2563 293 220 --3012 346 274 --2792 328 253 --2766 328 253 --2753 321 252 --3515 333 259 ---


31 (MPa) 2.51 2.72 1.5 2.1 --2.12 1.5 1.8 --2.10 1.5 1.6 --1.31 2.14 1.3 1.1 --1.60 1.2 1.5 --1.41 1.1 1.1 --1.21 1.1 1.0 --1.63 1.2 1.1 --0.80 1.2 1.3 ---

04 (K) 8390 8532 833 866 885 8221 801 852 849 8152 763 837 813 8202 7897 704 714 730 8341 848 893 931 8463 865 924 989 8321 866 921 989 8712 857 930 995 8821 862 947 1017

Page 35 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60



Figure Captions Figure 1. Molecules examined in this study. MW denotes the molecular weight of the molecule. Figure 2. Optimized molecular structures for C30H18 and C30H18NS. Gray, white, blue, and yellow atoms represent carbon, hydrogen, nitrogen, and sulfur, respectively. Figure 3. Particle number probability distributions from grand canonical Monte Carlo simulation for C30H26 at saturation conditions. Figure 4. Schematic of the process used to trace saturation curves. Open circles represent the results from direct grand canonical simulations. The dashed and solid curves represent the guess and converged, respectively, relationships associated with the isothermal-isobaric temperature expanded ensemble approach. Figure 5. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C30H18, C30H26, C30H34, and C30H48. Solid and dashed curves represent simulation data obtained with GAFF and CGenFF, respectively. Diamond, circle, and × symbols correspond to critical points for GAFF, CGenFF, and the NMRR correlation, respectively. Up and down triangles correspond to vapor pressures and enthalpies of vaporization, respectively, generated with the NMRR correlation. Squares represent predictions from the Reid et al. correlation. Figure 6. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C30H18, C34H26, C39H36 (S), and C39H36 (L). Curves and symbols are defined in the same manner as in Figure 5. Figure 7. Temperature dependence of saturated densities (top panel), vapor pressure (middle panel), and enthalpy of vaporization (bottom panel) for C39H36 (L), C38H35N, and C37H33NS. Curves and symbols are defined in the same manner as in Figure 5.



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 45

Liquid-vapor phase behavior of asphaltene-like molecules Figure 8. Normal boiling temperature. For a given molecule, the bars from left to right represent GAFF (solid blue), CGenFF (red diagonal hatch, not present for all molecules), MG (green square hatch), NMRR (violet horizontal hatch), aSB (maroon checkerboard).


Page 37 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 1 254x190mm (300 x 300 DPI)



1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Figure 2 254x190mm (300 x 300 DPI)


Page 38 of 45

Page 39 of 45

ln Π(N)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


0

-5

-10

950 K -15

-20

900 K

-25 0

20

40

60

80

N

100


120

140


Psat (MPa)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

10

10

10

0

Page 40 of 45

TEE - final TEE - guess GC - direct

-1

-2

-3

10 1.0

1.1

1.2

1.3

1.4

-1

1000 / T (K ) ACS Paragon Plus Environment

1.5

Page 41 of 45


1200

T (K)

1000 800 600 400

0

C30H18

C30H26

C30H34

C30H48

200

400

600

800

3

1000

1200

ρ (kg/m ) Psat (MPa)

10 10 10 10

0

-2

-4

-6

0.8

1.0

1.2

1.4

1.6

-1

1.8

2.0

1000/T (K ) ΔHvap (kJ/mol)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

150

100

50 500

600

700

800

900

Paragon Plus Environment TACS (K)

1000


1200

T (K)

1000 800

C30H18 C39H36 (L)

600 400

0

200

C34H26 C39H36 (S) 400

600

800

3

1000

1200

ρ (kg/m ) Psat (MPa)

10 10 10 10

0

-2

-4

-6

0.8

1.0

1.2

1.4

1.6

-1

1.8

2.0


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

150

100

50 500

600

700

800

900


1000

Page 42 of 45

Page 43 of 45


1200

T (K)

1000 800

C39H36 (L)

600

C37H33NS

400

Psat (MPa)

10 10 10 10

0

200

C38H35N

400

600

3

ρ (kg/m )

800

1000

1

0

-1

-2

10 10

-3

-4

10

-5

0.8

1.0

1.2

-1

1.4

1.6


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

140 120 100 80 60

600

700

800

900


1000


900

Aromaticity

800 700

Normal boiling point , Tb (K)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

C30H18

1000

C30H26

C30H34

C30H48

C39H36 (L)

C39H36 (S)

Alkyl chains

900 800 700 C30H18

C34H26

1100 Heteroatoms 1000 900 800 C39H36 (L)

C38H35N

C37H33NS


Page 44 of 45

Page 45 of 45

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TOC Graphic 254x190mm (96 x 96 DPI)


LiquidâVapor Phase Behavior of Asphaltene-like Molecules

Recommend Documents

LiquidâVapor Phase Behavior of Asphaltene-like Molecules