Vapor Liquid Equilibria of Hydrofluorocarbons Using Dispersion

Subscriber access provided by UNIV OF NEBRASKA - LINCOLN

Article

Vapor Liquid Equilibria of Hydrofluorocarbons using Dispersion-Corrected and Nonlocal Density Functionals Himanshu Goel, Charles L Butler, Zachary W Windom, and Neeraj Rai J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.6b00305 • Publication Date (Web): 13 Jun 2016 Downloaded from http://pubs.acs.org on June 14, 2016

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.

Journal of Chemical Theory and Computation 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 33

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

Journal of Chemical Theory and Computation

Vapor Liquid Equilibria of Hydrofluorocarbons using Dispersion-Corrected and Nonlocal Density Functionals Himanshu Goel,† Charles L. Butler,†,‡,¶ Zachary W. Windom,†,¶ and Neeraj Rai∗,† †Dave C. Swalm School of Chemical Engineering, and Center for Advanced Vehicular Systems, Mississippi State University, Mississippi State, MS, USA ‡East Mississippi Community College, Scooba, MS, USA ¶These authors contributed equally E-mail: [email protected] Phone: +1 (662) 3250790. Fax: +1 (662) 3252482

1

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

Abstract Recent developments in dispersion corrected and nonlocal density functionals are aimed at accurately capturing dispersion interactions, a key shortcoming of local and semilocal approximations of density functional theory. These functionals have shown significant promise for dimers and small clusters of molecules as well as crystalline materials. However, their efficacy for predicting vapor liquid equilibria is largely unexplored. In this work, we examine the accuracy of dispersion-corrected and nonlocal van der Waals functionals by computing the vapor liquid coexistence curves (VLCCs) of hydrofluoromethanes. Our results indicate that PBE-D3 functional performs significantly better in predicting saturated liquid densities than the rVV10 functional. With PBE-D3 functional, we also find that as the number of fluorine atoms increase in the molecule, the accuracy of saturated liquid density prediction improves as well. All the functionals significantly under predict the saturated vapor densities, which also result in an under prediction of saturated vapor pressure of all compounds. Despite the differences in the bulk liquid densities, the local microstructure of the liquid CFH3 and CF2 H2 are relatively insensitive to the density functional employed. For CF3 H, however, rVV10 predicts slightly more structured liquid than PBE-D3 functional.

2


Page 2 of 33

Page 3 of 33

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


1

Introduction

Although relatively weak when compared to covalent and hydrogen bonds, van der Waals (vdW) interactions play an important role in chemical and biological systems as they govern thermophysical properties of molecular systems. 1–4 The weak attractive interactions (also referred to as London dispersion interactions 5 ) arise from dipole-dipole (and higher order) polarizabilities of molecules. 6,7 Empirical force fields such as CHARMM, 8 AMBER, 9 OPLS, 10 and TraPPE 11 are currently the preferred choice for simulating condensed phase systems that are governed by dispersion interactions. However, advances in computational hardware and numerical algorithms allow the use of electronic structure calculations, in particular, density functional theory (DFT), to model the potential energy of condensed phase systems. Despite growing popularity, there are impediments in using DFT to model such systems as the commonly employed exchange-correlation functionals based on local density approximation (LDA), 12,13 semilocal (GGA), 14–17 and hybrid functionals, 18,19 lack long range electron correlation necessary to account for dispersion interactions at medium to large intermolecular distances. 20–24 There has been significant effort in the field to remedy this shortcoming either by purely empirical means 25–29 or by methods based on rigorous theoretical basis. 30–35 Among the empirically derived functionals, the “dispersion-corrected” DFT-D2 25 and DFT-D3 26 methods developed by Grimme and co-workers are quite successful and widely used. The main advantage of these methods is the negligible additional computational cost as it adds two-body or three-body dispersion terms to the standard Kohn Sham (KS)-DFT energy. At the other end of the spectrum (cost and theoretical foundation) are the nonlocal vdW density functionals. 30,31,36,37 The nonlocal correlation functional methods explicitly add a nonlocal term to the local or semilocal functionals. In this work, we have used the recent rVV10 functional, 33 which is a revised version of VV10 nonlocal density functional by Vydrov and Van Voorhis. 32 The VV10 functional is particularly successful and accurate in describing weakly interacting systems. 32,33 3



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

Although there have been numerous studies 20,35,38–41 that have looked at the efficacy of DFT functionals for describing vdW interactions present in small molecular clusters or crystalline materials, their efficacy for predicting vapor liquid equilibria (VLE) is largely unexplored. Understanding and predicting VLE of small organic molecules from first principles can have a profound impact on separation process such as sorption of organic molecules in nanoporous materials. In a significant contribution, McGrath and coworkers 42 carried out first principles Monte Carlo (FPMC) simulations of water with BLYP functional 15,17 and triple zeta basis set. They found that the saturated liquid density, the normal boiling point, and the critical temperature were significantly underpredicted. In 2011, McGrath et al. 43 used DFT-D2 25 with BLYP functional to study VLE of methane and methanol with FPMC simulations. In this work, as opposed to water simulations, BLYP with D2 correction significantly overestimated saturated liquid density and critical temperature for both molecules. Based on these simulations, it appears that the BLYP functional, without dispersion correction, it does not adequately account for dispersion interactions, whereas with D2 corrections overestimates the strength of these interactions in the condensed phase. Since the work of McGrath et al., 43 Grimme’s updated DFT-D3 26 and rVV10 33 functionals have been quite successful in describing dispersion interactions. Thus, it is plausible that DFT-D3 and rVV10 non local functional can predict vapor liquid equilibria with greater accuracy. The aim of the present work is to assess the accuracy of KS-DFT 44 by using different models that can account for dispersion interactions (PBE-D3, BLYP-D3, and rVV10) for predicting the vapor liquid coexistence curves (VLCC) of hydrofluoromethanes via FPMC simulations. We selected hydrofluoromethanes as hydrofluorocarbons(HFCs) are extensively used as a refrigerant in refrigerators and air conditioning units, fire suppression agent, propellant in the metered dose inhalers, etching agent, and blowing agent purposes. 45,46 The atmospheric concentration of HFCs is rapidly increasing and causing great environmental concern due to their greenhouse potential. 45,47 Due to their technological and environmental relevance, HFCs have been the focus of a number of experimental 48–51 and simulation stud-

4


Page 4 of 33

Page 5 of 33

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: Molecular structure of CF3 H (A), CF2 H2 (B), and CFH3 (C). The molecular frame of reference used for spatial distribution functions is also shown. The black, white and cyan colors represent carbon, hydrogen and fluorine atoms, respectively.

ies. 52–62 Furthermore, previous systems considered for FPMC had hydrogen bonding (water and methanol) and dispersion interactions (methane) as the dominant interactions. 42,43 HFCs selected here belong to a different class as these molecules, in addition to dispersion interactions, will have significant permanent dipole-dipole interactions, where upon thermal averaging, the leading order term goes as r−6 . 63 In the next section, simulation details for the first principles Monte Carlo simulations for hydrofluoromethanes are described. In Section 3, dipole moments, dimers energetics, structural analysis of the liquid phase, VLCCs, Clausius-Clapeyron plots, and critical properties are discussed, and Section 4 presents concluding remarks.

2

Simulation Details

The canonical version of Gibbs ensemble 64 Monte Carlo (GEMC) simulations were carried out to obtain VLCC for hydrofluoromethanes (CF3 H, CF2 H2 , and CFH3 ). Figure 1 shows the molecular structure of all the molecules considered in this study. GEMC simulations employ two separate periodic simulation boxes to represent vapor and liquid phases that are connected thermodynamically via an unified partition function. 65 This arrangement allows determination of VLE quantities such as saturated liquid density, saturated vapor density, and enthalpy of vaporization in a single simulation. The configurational part of phase space is efficiently sampled by translational, conformation change, rotational, volume, and swap

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

(particle exchange) moves. The translation, rotation and conformational change moves enable thermal equilibration, coupled volume changes (to maintain overall volume constant) enable mechanical equilibrium, and particle exchange between two boxes equilibrates chemical potentials. The volume and swap moves also provide the thermodynamic connection between the two separation periodic boxes. Grimme’s updated method DFT-D3 26 with PBE 16 and rVV10 non-local van der Waals density functionals 33 were used to simulate VLCCs for all of the hydrofluoromethanes. For fluoromethane, DFT-D3 was also examined with Becke-Lee-Yang-Parr (BLYP) 15,17 exchange and correlation functional. For rVV10 functional simulations, the VV10 van der Waals density funtional was added to the refitted version of Perdew-Wang 86 (PW86R), 36 which was used as the exchange part. The Perdew, Burke, and Ernzerhof (PBE) 16 generalized gradient approximation was used as the correlation part. All simulations were carried out using the Monte Carlo suite of CP2K simulation package. 66 The Quickstep module 67 was used to compute the interaction energies from first principles using Kohn-Sham density functional theory (KS-DFT). In all simulations, TZV2P basis set 67 and the pseudo potentials of Godecker et al. 68,69 with an auxiliary basis plane wave cutoff of 600 Ry were used. The PW grid was used to compute the derivatives for exchange-correlation. Given extremely high computational cost of GEMC simulations, approximate bias potentials were employed to increase the sampling efficiency. 70–73 We have used 16 moves with the approximate potential, before using DFT energies to accept or reject the move. Details of generating bias potentials and potential parameters are provided in the Supporting Information. The initial configurations for all hydrofluoromethanes were prepared using the Minnesota version of Monte Carlo for Complex Chemical Systems (MCCCS-MN) code. 74 The system consisted of 48 molecules. At each temperature of interest, the density was set to the experimental saturated liquid density and the volume of the vapor box was set identical to the liquid box containing a single molecule. 42 The reference cell for all simulations was set to 10% larger than the cell length of the box. A fixed reference cell length is necessary due

6


Page 6 of 33

Page 7 of 33

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: Configurations of hydrofluoromethane dimers considered for dimer energetics.

to the fluctuating box volume in GEMC simulations. 42 All of the simulations were run for approximately 600-800 Monte Carlo cycles, where each cycle consists of 48 moves. The initial 300-400 cycles were used for equilibration and the remaining cycles were considered to be the production run. The production run cycles were divided into blocks of 75 cycles each to compute averages and standard deviations. Each run of 600-800 MC cycle consisted of only one independent run. The 600-800 MC cycles were accomplished by restarting the simulation with the coordinates from the previous run. Typically, we were able to run 50-100 cycles in 48 hours with 80 processor cores (Xeon E5-2680 v2 at 2.80 GHz, 20 cores on each node, and nodes are connected by FDR InfiniBand interconnect). In total, 35 sets of simulations (each with 600-800 MC cycles) were carried out to explore the phase equilibria for all compounds. The density scaling law (with critical exponent= 0.325) and the law of rectilinear diameter were used to calculate critical temperature and critical density for each compound from VLCCs. 75–77 Also, the normal boiling point for each compound was calculated via ClausiusClapeyron equation and three or four temperatures close to the normal boiling point were considered. The liquid structure was analyzed using radial distribution function, angular radial distribution function (ARDF), and three-dimensional spatial distribution function (SDF). The details of ARDF and SDF calculations can be found in our previous work. 78–80

7



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 8 of 33

Table 1: Dipole Moments in Debye Units for Hydrofluoromethanes with Different DFT Functionals Functional/Basis PBE/PW BLYP/PW PBE-D3/PW BLYP-D3/PW rVV10/PW PBE/avTZ BLYP/avTZ PBE-D3/avTZ BLYP-D3/avTZ Experimental. 94 PW: Plane wave basis avTZ: aug-cc-pVTZ

CFH3 1.891 2.009 1.893 2.009 1.970 1.767 1.866 1.767 1.873 1.850

% Error 2.2 8.6 2.4 8.6 6.5 -4.5 0.9 -4.5 1.3 -

CF2 H2 1.941 2.049 1.943 2.054 2.018 1.854 1.948 1.857 1.954 1.970

% Error -1.5 4.0 -1.4 4.3 2.5 -5.9 -1.1 -5.8 -0.9 -

CF3 H 1.592 1.661 1.594 1.665 1.65 1.543 1.608 1.545 1.613 1.650

% Error 3.5 0.7 -3.4 0.9 0.0 -6.5 -2.6 -6.4 2.3 -

We also compared the DFT potential energy curves (PECs) for dimers with PECs obtained at coupled cluster single, double, and perturbative triple excitations (CCSD(T)) 81–86 level of theory. The Dunning’s augmented correlation-consistent polarized valence triple zeta basis set of contracted Gaussian functions (aug-cc-pVTZ) 87,88 with second order MP2 theory 89,90 is used to determine the minimum energy orientation of dimers, and to optimize the monomer geometries. The configurations used for different molecules are shown in Figure 2. The geometries were kept rigid during the potential energy scans. The intermolecular binding energies were corrected for basis set superposition error (BSSE) using the counterpoise approach. 91 We used Gaussian 92 and Molpro 93 software packages for CCSD(T) and other dimer calculations.

3 3.1

Results and Discussion Dipole Moments and Dimer Energetics

As all the molecules considered here have significant dipole moments, permanent dipoledipole interactions are expected to play an important role. The gas phase dipole moments for CF3 H, CF2 H2 and CFH3 with different functionals are reported in Table 1. Experimentally, the hydrofluoromethanes exhibit dipole moments in the following order CF2 H2 >CFH3 8


Page 9 of 33

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


˚), Interaction Energy (E int , kJ/mol), and Basis Set Table 2: Location of Well Depth (R 0 , A Superposition Error (BSSE, kJ/mol) of Hydrofluoromethane at the minima in the Dimer Potential Energy Curve for DFT Functionals Compared to the CCSD(T)/aug-cc-pVTZ Level of Theory Configuration 1 Functional PBE-D3/PW BLYP-D3/PW PBE-D3/avTZ BLYP-D3/avTZ rVV10 CCSD(T) Configuration 2 Functional PBE-D3/PW BLYP-D3/PW PBE-D3/avTZ BLYP-D3/avTZ rVV10 CCSD(T)

R0 4.45 4.375 4.40 4.30 4.375 4.375 R0 2.50 2.45 2.50 2.45 2.40 2.47

CFH3 -E int 6.705 6.212 6.485 6.192 6.295 6.621 CFH3 -E int 10.738 10.866 10.460 10.460 10.329 9.870

BSSE 0.715 1.30 0.230 0.272 0.773 1.037

R0 4.0 3.875 4.0 3.875 4.0 3.95

BSSE 1.027 1.73 0.392 0.436 1.233 1.922

R0 2.70 2.55 2.70 2.60 2.60 2.65

CF2 H2 -E int 7.864 7.135 7.531 6.778 7.072 7.230 CF2 H2 -E int 11.928 12.423 11.715 11.966 11.637 12.712

BSSE 1.340 0.036 0.338 0.388 0.783 1.780

R0 4.15 4.0 4.15 4.05 4.10 4.15

BSSE 0.286 1.278 0.641 0.718 1.284 2.780

R0 2.82 2.71 2.8 2.7 2.70 2.70

CF3 H -E int 5.478 4.691 5.230 5.188 4.127 4.260 CF3 H -E int 9.018 9.208 8.911 8.953 8.838 9.480

BSSE 1.287 1.781 0.602 0.647 0.767 2.420 BSSE 1.788 1.938 0.652 0.736 0.363 3.123

>CF3 H. 94 All functionals considered here predict the correct trend in molecular dipole moment. The dipole moment values for underlying functionals (PBE, BLYP) with and without dispersion corrections are very close to each other. Adding dispersion correction terms to these underlying functionals helps in capturing the long range attractive forces, but they do not significantly alter gas phase dipole moments. In comparison between plane wave (PW) and avTZ (aug-cc-pVTZ) basis set for BLYP or PBE functional, the dipole moment values are slightly larger for PW basis set. With both functionals (PBE/BLYP), the difference is around 3%, 5% and 7% for CF3 H, CF2 H2 and CFH3 , respectively. When using PW basis set, the PBE-D3 functional performs better than BLYP-D3 and rVV10 functionals. For all three compounds, the mean unsigned percentage error for PBE-D3 functional is 2.3% as compared to 4.5% and 3% with BLYP-D3 and rVV10 functionals, respectively. For CF2 H2 and CF3 H molecules, most functionals considered here predict reasonable dipole moments as compared with experimental values, whereas the error is higher for CFH3 molecule. As the two–body terms contribute most to the n–body system potential energy, we

9



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

examine dimer potential energy curves. For the energy calculation, two sets of different hydrofluoromethane configurations were considered. We performed exploratory calculations to identify the most stable dimer configuration that maximizes dipole-dipole interactions for each of these molecules at the MP2 level of theory. These configurations (CONF1) are shown in the Figure 2 where r is the distance between carbon atoms. In addition, the dimer configurations with direct H-F interactions were also considered in configuration 2, where r is the distance between hydrogen and fluorine atom (Figure 2). We have considered the dimer configurations for direct H-F interactions from Kryachko and Scheiner’s work 62 where they have examined the hydrogen bonding between different pairs of hydrofluoromethanes at MP2 level of theory. For configuration 1, Figure 3 shows the potential energy curves (PECs) for each hydrofluoromethane dimer configuration with CCSD(T)/avTZ level of theory, 87 PBE-D3/PW, rVV10/PW, BLYP-D3/PW, BLYP-D3/avTZ, and PBE-D3/avTZ functionals. The equilibrium bond lengths, well depths, and basis set superposition errors (BSSEs) for both configurations are presented in Table 2. The binding energy obtained from different density functionals lie within 1 kJ/mol of CCSD(T) calculations for both configurations of hydrofluoromethanes. The well depth and equilibrium bond length obtained from plane wave and avTZ basis set with GGA functional are in good agreement with each other. In general, the well depth and equilibrium bond length obtained from rVV10 functional matche closely with the results obtained at the CCSD(T) level of theory. Similar results were obtained for rVV10 density functional for dispersion-bound, mixed, and hydrogen-bonded complexes (S22 set), 33 rare gas dimers and solids. 95 On the other hand, the PBE-D3 functional slightly overestimates the equilibrium bond length for most of the hydrofluoromethanes.

10


Page 10 of 33

Page 11 of 33

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 3: Potential energy curves of hydrofluoromethanes with different functionals for configuration 1. Figure A, B, and C present the potential energy curves for CF3 H, CF2 H2 and CFH3 , respectively. The black star, green circle, red cross, magenta plus, blue diamond, and brown triangle symbols represent data obtained from the CCSD(T)/avTZ level of theory, PBE-D3/PW, rVV10/PW, BLYP-D3/PW, BLYP-D3/avTZ, and PBE-D3/avTZ functionals, respectively. The dashed lines of the corresponding color are guide to the eye. 11



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 4: Liquid phase radial distribution functions for all hydrofluoromethanes. The solid, dashed, and dashed dotted lines represent plots for CF3 H, CF2 H2 and CFH3 , respectively. The red, black, and green colors depicts PBE-D3, rVV10 and BLYP-D3 functionals, respectively. The RDF plots are plotted at T = 265 K for CF3 H, T = 315 K for CF2 H2 and T = 270 K for CFH3 . The offset for CF2 H2 and CFH3 is 1 and 2 units in the y direction, respectively.

3.2

Liquid Structure

We employed radial distribution function (RDF) to investigate the effect of the functionals on the liquid microstructure. The radial distribution functions (RDFs) for all hydrofluoromethanes (between carbon-carbon, hydrogen-hydrogen, fluorine-fluorine, and fluorinehydrogen) with different functionals are shown in Figure 4. As the experimental data is rather limited for these molecules, we compared the liquid for CF2 H2 and CF3 H with previous atomistic simulation based on pair-wise additive potential. 52 In this work, the RDF plots were computed at 270 K, 315 K, and 265 K for CFH3 , CF2 H2 and CF3 H, respectively.

12


Page 12 of 33

Page 13 of 33

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 5: Potential H–F hydrogen bonding of two CF2 H2 molecules which results energetically highly favorable configuration.

L´ısal and Vacek’s 52 RDF plots for CF2 H2 and CF3 H were determined at similar temperatures (313.15 K, and 266.6 K for CF2 H2 and CF3 H, respectively). Additionally, we have also computed center of mass-center of mass radial distribution functions for CF2 H2 and CF3 H at 315 K and 265 K, respectively. The location of the first coordination shell peak are at about 4.46 ˚ A with PBE-D3 and 4.23 ˚ A with rVV10 for CF2 H2 and for CF3 H, the location is at 4.68 ˚ A with PBE-D3 and 4.27 ˚ A with rVV10. L´ısal and Vacek 52 have determined center-center radial distribution functions for CF2 H2 and CF3 H. The location of the first coordination peaks are at about 4.4 ˚ A for CF2 H2 and 4.7 ˚ A for CF3 H. The PBE-D3 results for CF2 H2 and CF3 H are in very close agreement with Lisal and Vacek’s 52 work. The C-C RDFs first coordination shell peaks are at about 4.30 ˚ A (with PBE-D3), 4.29 ˚ A (with rVV10) for CF2 H2 , and 4.6 ˚ A (with PBE-D3), 4.22 ˚ A (with rVV10) for CF3 H. The liquid microstructure obtained with the first principles calculations with PBE-D3 is in good agreement with Lisal and Vacek’s 52 work despite significant over prediction of the saturated liquid densities (vide infra). The C–C and F–H RDF plots for CF2 H2 with rVV10 show a small peak located at 3 and 1.7 ˚ A, respectively. This feature is not observed for PBE-D3 functional or previous simulations based on molecular mechanics force fields. As the first prominent peak in different RDFs of CF2 H2 suggest the interlocking mechanism proposed initially by Lowden and Chander 96 (for liquid carbon tetrachloride) and observed by L´ısal and Vacek for CF2 H2 , we looked 13



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 6: Angular radial distribution functions (ARDF) for CF3 H (A), CF2 H2 (B) and CFH3 (C).

at the simulation trajectory to identify the configuration that corresponds to this anomalous peak. One such configuration is shown in Figure 5. This configuration does not conform to the interlocking molecular arrangement where one of the fluorine atoms sits in one of the four cavities of the neighboring molecules due to the near tetrahedral geometry. 52,96 It appears that rVV10 functional enhances F–H interactions that result in stronger F–H intermolecular hydrogen bonds. This could also be the reason for over prediction of the saturated liquid densities (vide infra). Further important structural details about positioning of the molecules can be obtained from angular radial distribution functions (ARDFs). Since the ARDF resolves angles in

14


Page 14 of 33

Page 15 of 33

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


addition to the radial distance (as in RDFs), they can provide information about the relative orientation of the molecules in the liquid phase. Figure 6 shows ARDF for CFH3 , CF2 H2 and CF3 H molecules (with PBE-D3 functional) at 270, 315, and 265 K respectively. For the ARDF, distance is computed between the carbon atoms for all hydrofluoromethanes. The angles were defined based on the vectors that closely align with the dipole moment vector of the molecule. Thus, for CF3 H, the angle is calculated between vectors from hydrogen to carbon, and carbon to fluorine for CFH3 . For CF2 H2 , the angle is defined between the bisector of the FCF angle. From the ARDF plot, all the hydrofluoromethane molecules show a strong inclination for parallel or antiparallel orientation in the first solvation shell. The fluoromethane ARDF plot shows considerably higher peaks for antiparallel orientation as compare with other molecules. The fluoromethane shows one prominent peak at r ≈ 3.41 ˚ A and θ ≈ 178°, and another one at r ≈ 3.81 ˚ A and θ ≈ 1.25°, respectively. This denotes that a fluoromethane molecule desires to be in antiparallel orientation such that fluorine can hydrogen bond with hydrogen atoms. A very similar pattern is observed with trifluoromethane molecules. There are many tall peaks around r ≈ 4.46 ˚ A and θ ≈ 179 ° and a few of them at r ≈ 4.04 ˚ A and θ ≈ 1.25 ° as well, which again signifies parallel or antiparallel orientation and rearrangement of fluorine and hydrogen atoms due to hydrogen bonding. In difluoromethane, one would expect to have a large extent of hydrogen bonding due to a higher dipole moment. There are two fluorine atoms on one side of carbon and two hydrogen atoms on the other side of carbon. The Figure 6 B shows that the CF2 H2 molecule show somewhat smaller peaks in the first solvation shell as compared with the other two molecules. There are a lot of smaller height peaks present in the range of 2° to 30° around 4.1 ˚ A. Figure 7 shows spatial distribution functions(SDFs) for liquid phase of CFH3 , CF2 H2 and CF3 H molecules (with PBE-DE functional) at 270, 315, and 265 K, respectively. The spatial distribution function is defined as three-dimensional density distribution of molecules around a central molecule. The surfaces shown in the figure correspond to the 25% most

15



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 7: Spatial distribution functions (SDF) for the liquid phase of CF3 H (A), CF2 H2 (B) and CFH3 (C).The 25 % most likely positions around the particular molecule is shown.

likely position of molecules in the first solvation shell. Figure Figure 7 A and Figure 7 C shows higher density distribution around the single hydrogen and single fluorine atom, respectively. This is due to the hydrogen bonding, where single hydrogen or fluorine atoms are encircled by three fluorine or hydrogen atoms. Figure 7 B indicates that the most preferred location of CF2 H2 molecules is along the F–C–F and H–C–H angle bisectors. The spatial distribution function also highlight the fact that stronger F–H interactions can be achieved in the first solvation shell when molecules are off-axis (molecular dipole vector) in the case of CFH3 and CF3 H.

3.3

Vapor Liquid Co-existence Curves

The GEMC simulations were carried out with three different functionals (BLYP-D3, PBED3 and rVV10) for CFH3 . Figure 8 shows the computed vapor-liquid coexistence curves (VLCCs) and Clausius-Clapeyron plots for CFH3 with different functionals. PBE-D3 and BLYP-D3 predictions for saturated liquid density are more accurate than rVV10. The PBED3 functional performs reasonably well across the range of T = 270-300 K as compared to BLYP-D3. The density overestimation for this range is approximately 20%. Overall, the average liquid densities obtained from PBE-D3 are overestimated by 29%, and with the BLYP-D3 and rVV10 functionals, densities are overestimated by 32% and 38%, respectively. 16


Page 16 of 33

Page 17 of 33

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 8: Vapor-liquid coexistence curves and Clausius-Clapeyron plots for CFH3 . The solid green lines depict the experimental data 97 and stars represent the experimental critical point for CFH3 . The black circles, red triangles, and blue squares represents GEMC simulation data with PBE-D3, BLYP-D3, and rVV10 functionals, respectively. The filled symbols represent the critical temperature for the corresponding functional. The dashed lines shows linear fits of the simulated data for the corresponding functional based on the color.

Clearly, rVV10 functional significantly overestimates the strength of intermolecular interactions in the condensed phase. McGrath et al. 43 also observed consistent error greater than 30% in the saturated liquid densities for methane and methanol with BLYP-D3 functionals. For the vapor phase, on the other hand, all the functionals significantly under predict the vapor density. Furthermore, all the functionals incorrectly predict the coefficient of thermal expansion due to the incorrect curvature of liquid and the vapor coexistence lines. Given very low vapor densities at all subcritical temperatures, we have used the ideal gas law to estimate saturated vapor pressure from the corresponding densities obtained from VLCCs. In Clausius-Clapeyron plot, the saturated vapor pressures fall significantly below the experimental values. 97 The PBE-D3 functional yields saturated vapor pressures below the other two functionals. As the performance of PBE-D3 and BLYP-D3 is similar for the saturated liquid part of the VLCC and extremely high computational cost of each VLCC, we decided to carry 17



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 9: Vapor-liquid coexistence curves and Clausius-Clapeyron plots for CF2 H2 . The solid orange lines depict the experimental data 97 and stars represent the experimental critical point for CF2 H2 . The black diamonds and blue triangles represent GEMC simulation data with PBE-D3 and rVV10 functional. The filled symbols represent the critical temperature for the corresponding functional. The dashed lines shows linear fits of the simulated data for the corresponding functional based on the color.

out simulations for the rest of the molecules with just PBE-D3 and rVV10. The computed vapor-liquid coexistence curves and Clausius-Clapeyron plots for CF2 H2 and CF3 H with PBE-D3 and rVV10 functional is shown in Figure 9 and 10. For CF2 H2 , both functional predictions are very close to each other and overestimate the saturated liquid density by approximately 25% at higher temperatures and 15% at lower temperatures. The percentage errors for saturated liquid density increase significantly at higher temperatures. In Figure 10, the saturated liquid density for CF3 H with PBE-D3 functional shows promising results when compared with experimental data. The mean unsigned percentage error of saturated liquid density with PBE-D3 functional is around 5%, whereas rVV10 overestimates by 28%. The liquid density with rVV10 is greatly overestimated for CF3 H and CFH3 when compared with the PBE functional. As with CFH3 , the saturated vapor pressures of CF3 H and CF2 H2 for both rVV10 and PBE-D3 fall considerably below experimental data. The PBE-D3 functional for CF3 H and CF2 H2 yields saturated vapor pressures higher than rVV10 functional 18


Page 18 of 33

Page 19 of 33

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 10: Vapor-liquid coexistence curves and Clausius-Clapeyron plots for CF3 H. The solid green lines depict the experimental data 97 and stars represent the experimental critical point for CF3 H . The black circles and blue squares represent GEMC simulation data with PBE-D3 and rVV10 functional. The filled symbols represent the critical temperature for the corresponding functional. The dashed lines shows linear fits of the simulated data for the corresponding functional based on the color.

simulations. When describing hydrofluoromethanes with the PBE-D3 functional, we note that as the number of fluorine atoms increases, the accuracy of the saturated liquid density and vapor pressure prediction improves as well. Furthermore, as the high dipole moments of compounds (vide infra) increases, PBE-D3 predicts a steeper saturated liquid line which results in larger error in the predicted critical temperature. The GGA functionals with dispersive correction DFT-D3, although significantly better than rVV10, are still not accurate enough to capture dispersive interactions satisfactorily for modeling vapor liquid equilibria. The primary purpose of simulating vapor liquid equilibria is to compute the critical point due to its fundamental role in developing equations of state for thermodynamic modeling. 98,99 For all compounds, the predicted critical temperatures (Tc ), critical density (ρc ) and normal boiling points (Tb ) are listed in Table 3. For CFH3 , the error in the critical temperature is around 5-6% with the BLYP-D3 and rVV10 functionals as compared to 20% with PBE-D3. The good agreement with the experimental critical points for BLYP-D3 and 19



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 33

Table 3: Critical Temperature (TC ), Critical Density (ρC ), and the Normal Boiling Point (TB ) for Hydrofluoromethanes Obtained From First Principles Monte Carlo and Classical Simulations. Experimental Data is Taken from the NIST Chemistry Webbook. 97 The Numbers in the Parenthesis are the Standard Deviation Compound CFH3

CF2 H2

CF3 H

Functional/ Force Field PBE-D3 rVV10 BLYP-D3 TraPPE OPLS Exp. PBE-D3 rVV10 TraPPE OPLS Exp. PBE-D3 rVV10 TraPPE OPLS Exp.

Tc (K)

ρc (g/cm3 )

Tb (K)

379(13) 334(4) 337(3) 254.5(1.2) 206.3(1.1) 317.28 440(15) 440(30) 255(2) 217.0(1.2) 351.25 298(4) 314 (3) 228.0(1.2) 212(2) 299.3

0.28(0.04) 0.30(0.03) 0.280(0.013) 0.310(0.005) 0.310(0.005) 0.316 0.40(0.04) 0.36(0.05) 0.412(0.007) 0.420(0.007) 0.424 0.42(0.04) 0.43(0.04) 0.507(0.010) 0.510(0.012) 0.526

270(10) 210(20) 210(14) 159(1) 128(1) 195.0 260(5) 240(20) 163(1) 137(1) 221.4 205(8) 160(52) 147(1) 136(1) 189.0

rVV10 is fortuitous as although the overall error for the liquid line of VLCC is larger for these functionals, the lines have smaller slope than the experimental curve, especially at the higher temperatures. The PBE-D3 liquid has a smaller coefficient of thermal expansion than the real fluid which manifests in a higher critical temperature. The critical density is only slightly under predicted for all three functionals. We have obtained normal boiling point for all compounds by fitting to the Clausius-Clapeyron equation. The normal boiling point for CFH3 leads to an error of around 6% with BLYP-D3 and rVV10, whereas it is nearly 26% with PBE-D3. Similarly, the predicted critical temperature for CF2 H2 with PBE-D3 and rVV10 functional is somewhat close to each other and results in an overestimation of about 24% and 25%. The normal boiling point for CF2 H2 compound results in an error of around 17% and 8% with PBE-D3 and rVV10 functional. Interestingly, the predicted critical temperature for CF3 H with PBE-D3 functionals are in excellent agreement with the exper-

20


Page 21 of 33

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


imental data. The percentage error is less than 1% with PBE-D3 whereas it is nearly 5% with the rVV10 functional. However, the critical density for CF3 H with both functionals was underestimated in a substantial way due to the underestimation of the vapor density. This trend is apparent in nearly all of the molecules and functionals presented here. The error in the predicted normal boiling point for CF3 H is around 8% and 17% with PBE-D3 and rVV10 functionals, respectively. For the lower dipole moment compound (CF3 H) VLCCs, critical temperature, and boiling temperature obtained from the PBE-D3 functional matche well with the experimental results. For the higher dipole moment compounds (CF2 H2 , CFH3 ), the VLCCs, critical properties, and boiling temperature obtained from PBE-D3 functional show significant over prediction. It is common to assume that potential energy of a N -particle system can be considered to be pair-wise additive and contribution of three–body, four–body, and n–body terms can be neglected. 100 However, this may not always be correct as induction terms have been shown to have many-body terms as the leading term. 100 Although we have not generated complete potential energy surface for the dimer interaction, the potential energy curves presented in Section 3.1 can be useful as they represent energetically most favorable interactions. Thus, these configurations are expected to contribute most to the configurational integral. For dimer energetics, overall rVV10 functional performs better than other DFT functionals considered in this study, however, vapor liquid coexistence curves deviate significantly from the experimental data. The rVV10 functional performs well for dimers of S22 set, 32,33 rare gas and solids, 95 but there is still room for improvement when considering condensed phase properties perhaps by reparameterizing b parameter as done by Björkman et al. for solids. 101,102 On the other hand, the PBE-D3 functional qualitatively fails to reproduce the dimer energetics for CF3 H, but works better for VLE. When examining the VLCC results for different functionals in the context of dimer energetics, it appears a good match (with CCSD(T)) for the dimers energetics by a functional may not necessarily guarantee good performance for condensed phase properties of the compounds. This may be partly because n-body terms

21



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

(where n > 2) can be significant for certain classes of molecules. Given apparently large errors in the saturated liquid densities and vapor pressures, we looked at the ability of common pair-wise additive force fields to predict VLCC for the molecules considered here. There are a number of studies, 52–54,56,103–106 where each hydrofluoromethane compound was parameterized against experimental data. Potter et al. 107 developed transferable potential model for CF3 H, CF2 H2 , and CF4 compounds. Their force field matches exceptionally well for CF3 H, when compared with experiment, but the model could not do a good job for other compounds. Clearly, if a force field is parameterized for a particular molecule against a set of experimental data, it will perform well. However, the real benefit of the FPMC approach is the case where experimental data is not available or is difficult to obtain. Thus, we picked two popular force fields (OPLS 108 and TraPPE 59 ) that are not specifically parameterized for hydrofluoromethanes for predicting the VLCCs and critical properties but are expected to have better transferability. The Lennard Jones parameters for TraPPE force field were taken from Rai et al. work 59 and partial atomic charges were obtained according to the TraPPE-EH approach 109 by performing electronic structure calculation (at B3P86/6-311G* level). In case with OPLS force field, the non bonded parameters for carbon and fluorine were taken from Watkins and Jorgensen’s work 108 while the hydrogen parameters are from OPLS-AA force field. 10 The bonded parameters considered for the classical simulations and the calculated GEMC simulation data are given in the Supporting Information. For all compounds, the VLCCs simulation data obtained from both force fields is substantially under predicting the saturated liquid and vapor densities, and critical properties. In Table 3, the critical temperatures for all three compounds with both force fields is listed. For all three compounds, the error for critical temperature is around 24% and 34% with TraPPE and OPLS force fields, respectively. It appears that PBE-D3 results for VLE of hydrofluoromethanes are significantly better than what one would obtain when using classical force fields without any parameterization.

22


Page 22 of 33

Page 23 of 33

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


4

Conclusions

In this work, we have presented first principles Monte Carlo calculations for CF3 H, CF2 H2 ,and CFH3 with different dispersion-corrected and nonlocal density functionals. Vapor liquid coexistence curves, saturated vapor pressures, critical temperatures, normal boiling points, dipole moments, dimer energetics and structural properties (RDFs, ARDFs, SDFs) of hydrofluoromethanes were determined. The rVV10 non local van der Waals density functional predictions were not very encouraging for predicting condensed phase properties of hydrofluoromethanes. It may be possible to improve condensed phase properties by reparameterizing the b parameter in rVV10 functional that controls the range of short range interactions. 101,102 We find that the PBE-D3 functional works particularly well for predicting saturated liquid density as the number of fluorine atoms increase in hydrofluoromethanes. The PBE-D3 functional shows promising results for the lower dipole moment compound CF3 H where predicted accuracy of the saturated liquid densities is comparable to highly parameterized empirical force fields. All the functionals with plane wave basis set predict the gas phase dipole moment with less than 5% error. The RDF analysis for all hydrofluoromethanes with different dispersion correction were in close agreement with previously published simulation studies based on empirical force fields. The ARDF and SDF analyses provide detailed microscopic view of first coordination shell where parallel and antiparallel configurations are preferred.

5

Acknowledgements

This work is funded by the National Science Foundation (NSF) under grant number CHE - 1265872. This work used the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation grant number OCI-1053575. XSEDE provided the computing support from the Texas Advanced Computing Center (TACC) at The University of Texas at Austin and University of Tennessee and Oak Ridge National Laboratory’s Joint Institute for Computational Sciences. This research also used 23



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

resources of the National Energy Research Scientific Computing Center, a DOE Office of Science User Facility supported by the Office of Science of the U.S. Department of Energy under contract number No. DE-AC02-05CH11231. We also acknowledge computational resources at High Performance Computing Collaboratory in Mississippi State University. We would like to thank Matthew McGrath and Joost VandeVondele for helpful discussions. In addition, we thank Jacek Jakowski and Justin Oelgoetz for their help with porting CP2K onto stampede.

Supporting Information Available The detailed description of the approximate potential parameters, numerical data figures 8, 9, and 10, numerical VLCC data for classical force fields, and sample input files are provided in the Supporting Information. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Brutschy, B.; Hobza, P. Chem. Rev. 2000, 100, 3861–3862. (2) Buckingham, A. D.; Fowler, P. W.; Hutson, J. M. Chem. Rev. 1988, 88, 963–988. (3) Chalasinski, G.; Gutowski, M. Chem. Rev. 1988, 88, 943–962. (4) Jeziorski, B.; Moszynski, R.; Szalewicz, K. Chem. Rev. 1994, 94, 1887–1930. (5) London, F. T. Faraday Soc. 1937, 33, 8b–26. (6) Maitland, G. C. Intermolecular Forces: Their Origin and Determination; Oxford University Press, 1981; pp 57–76. (7) Stone, A. The Theory of Intermolecular Forces; Oxford University Press, 2013; pp 124–140. 24


Page 24 of 33

Page 25 of 33

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


(8) Vanommeslaeghe, K.; Hatcher, E.; Acharya, C.; Kundu, S.; Zhong, S.; Shim, J.; Darian, E.; Guvench, O.; Lopes, P.; Vorobyov, I.; Mackerell Jr., A. D. J. Comput. Chem. 2010, 31, 671–690. (9) Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2004, 25, 1157–1174. (10) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. J. Am. Chem. Soc. 1996, 118, 11225–11236. (11) Martin, M. G.; Siepmann, J. I. J. Phys. Chem. B 1998, 102, 2569–2577. (12) Perdew, J. P.; Wang, Y. Phys. Rev. B 1992, 45, 13244–13249. (13) Vosko, S.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200–1211. (14) Perdew, J. P. Phys. Rev. B 1986, 33, 8822–8824. (15) Becke, A. D. Phys. Rev. A 1988, 38, 3098. (16) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865. (17) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B 1988, 37, 785. (18) Handy, N. C.; Cohen, A. J. Mol. Phys. 2001, 99, 403–412. (19) Schmider, H. L.; Becke, A. D. J. Chem. Phys. 1998, 108, 9624–9631. (20) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2007, 3, 289–300. ˇ y, J.; Hobza, P. Phys. Chem. Chem. Phys. 2007, 9, 5291–5303. (21) Cern` (22) Arey, J. S.; Aeberhard, P. C.; Lin, I.-C.; Rothlisberger, U. J. Phys. Chem. B 2009, 113, 4726–4732. (23) Klimeˇs, J.; Michaelides, A. J. Chem. Phys. 2012, 137, 120901. 25



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

(24) Kristyán, S.; Pulay, P. Chem. Phys. Lett. 1994, 229, 175–180. (25) Grimme, S. J. Comput. Chem. 2006, 27, 1787–1799. (26) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. J. Chem. Phys. 2010, 132, 154104. (27) Zhao, Y.; Truhlar, D. G. Accounts Chem. Res. 2008, 41, 157–167. (28) Von Lilienfeld, O. A.; Tavernelli, I.; Rothlisberger, U.; Sebastiani, D. Phys. Rev. Lett. 2004, 93, 153004. (29) Aeberhard, P. C.; Arey, J. S.; Lin, I.-C.; Rothlisberger, U. J. Chem. Theory Comput. 2009, 5, 23–28. (30) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D. C.; Lundqvist, B. I. Phys. Rev. Lett. 2004, 92, 246401. ´ D.; Kong, L.; Lundqvist, B. I.; Langreth, D. C. Phys. Rev. B (31) Lee, K.; Murray, E. 2010, 82, 081101. (32) Vydrov, O. A.; Van Voorhis, T. J. Chem. Phys. 2010, 133, 244103. (33) Sabatini, R.; Gorni, T.; de Gironcoli, S. Phys. Rev. B 2013, 87, 041108. (34) Kannemann, F. O.; Becke, A. D. J. Chem. Theory Comput. 2009, 5, 719–727. (35) Kannemann, F. O.; Becke, A. D. J. Chem. Theory Comput. 2010, 6, 1081–1088. ´ D.; Lee, K.; Langreth, D. C. J. Chem. Theory Comput. 2009, 5, 2754–2762. (36) Murray, E. (37) Vydrov, O. A.; Van Voorhis, T. Phys. Rev. Lett. 2009, 103, 063004. (38) Rao, L.; Ke, H.; Fu, G.; Xu, X.; Yan, Y. J. Chem. Theory Comput. 2008, 5, 86–96. (39) Hohenstein, E. G.; Chill, S. T.; Sherrill, C. D. J. Chem. Theory Comput. 2008, 4, 1996–2000. 26


Page 26 of 33

Page 27 of 33

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


(40) Sedlak, R.; Janowski, T.; Pito´ nak, M.; Re´za´ć, J.; Pulay, P.; Hobza, P. J. Chem. Theory Comput. 2013, 9, 3364–3374. (41) Sorescu, D. C.; Byrd, E. F.; Rice, B. M.; Jordan, K. D. J. Chem. Theory Comput. 2014, 10, 4982–4994. (42) McGrath, M. J.; Siepmann, J. I.; Kuo, I.-F. W.; Mundy, C. J.; VandeVondele, J.; Hutter, J.; Mohamed, F.; Krack, M. J. Phys. Chem. A 2006, 110, 640–646. (43) McGrath, M. J.; Kuo, I.-F. W.; Ghogomu, J. N.; Mundy, C. J.; Siepmann, J. I. J. Phys. Chem. B 2011, 115, 11688–11692. (44) Kohn, W.; Sham, L. J. Phys. Rev. 1965, 140, A1133. (45) Tsai, W.-T. Chemosphere 2005, 61, 1539–1547. (46) Case, F. H.; Brennan, J.; Chaka, A.; Dobbs, K. D.; Friend, D. G.; Frurip, D.; Gordon, P. A.; Moore, J.; Mountain, R. D.; Olson, J.; Ross, R. B.; Schiller, M.; Shen, V. K. Fluid Phase Equilibr. 2007, 260, 153 – 163, 3rd Fluid Properties Challenge. (47) Montzka, S. A.; McFarland, M.; Andersen, S. O.; Miller, B. R.; Fahey, D. W.; Hall, B. D.; Hu, L.; Siso, C.; Elkins, J. W. J. Phys. Chem. A 2015, 119, 4439–4449. (48) Moghadasi, J.; Papari, M. M.; Mohammad-Aghaie, D.; Campo, A. B. Chem. Soc. Jpn. 2008, 81, 220–234. (49) Neuefeind, J.; Fischer, H.; Schröer, W. J. Phys-Condens Mat. 2000, 12, 8765. (50) Mort, K.; Johnson, K.; N Burgess, A. Mol. Phys. 2000, 98, 999–1003. (51) Santos, F.; Pai-Panandiker, R.; De Castro, C.; Mardolcar, U. IEEE Trans. Dielectr. Electr. Insul. 2006, 13, 503–511. (52) Lisal, M.; Vacek, V. Fluid Phase Equilibr. 1996, 118, 61–76.

27



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

(53) Böhm, H.; Meissner, C.; Ahlrichs, R. Mol. Phys. 1984, 53, 651–672. (54) Fermeglia, M.; Ferrone, M.; Pricl, S. Fluid Phase Equilibr. 2003, 210, 105–116. (55) Peguin, R. P.; Kamath, G.; Potoff, J. J.; da Rocha, S. R. J. Phys. Chem. B 2008, 113, 178–187. (56) Stoll, J.; Vrabec, J.; Hasse, H. J. Chem. Phys. 2003, 119, 11396–11407. (57) Oakley, M. T.; Do, H.; Hirst, J. D.; Wheatley, R. J. J. Chem. Phys. 2011, 134, 114518. (58) Paulechka, E.; Kroenlein, K.; Kazakov, A.; Frenkel, M. J. Phys. Chem. B 2012, 116, 14389–14397. (59) Rai, N.; Rafferty, J. L.; Maiti, A.; Siepmann, J. I. Fluid Phase Equilibr. 2007, 260, 199–211. (60) Caminati, W.; López, J. C.; Alonso, J. L.; Grabow, J.-U. Angew. Chem., Int. Ed. 2005, 117, 3908–3912. (61) Cabral, B. J. C.; Guedes, R. C.; Pai-Panandiker, R. S.; de Castro, C. A. N. Phys. Chem. Chem. Phys. 2001, 3, 4200–4207. (62) Kryachko, E.; Scheiner, S. J. Phys. Chem. A 2004, 108, 2527–2535. (63) Hirschfelder, J. O.; Curtiss, C. F.; Bird, R. B.; Mayer, M. G. Molecular theory of gases and liquids; Wiley New York, 1954; Vol. 26; pp 22–30. (64) Panagiotopoulos, A. Z. Mol. Phys. 1987, 61, 813–826. (65) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic press, 2001; Vol. 1; pp 201–223. (66) Hutter, J.; Iannuzzi, M.; Schiffmann, F.; VandeVondele, J. WIREs Comp. Mol. Sci. 2014, 4, 15–25. 28


Page 28 of 33

Page 29 of 33

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


(67) VandeVondele, J.; Krack, M.; Mohamed, F.; Parrinello, M.; Chassaing, T.; Hutter, J. Comput. Phys. Commun. 2005, 167, 103–128. (68) Goedecker, S.; Teter, M.; Hutter, J. Phys. Rev. B 1996, 54, 1703. (69) Hartwigsen, C.; Gœdecker, S.; Hutter, J. Phys. Rev. B 1998, 58, 3641. (70) Iftimie, R.; Salahub, D.; Wei, D.; Schofield, J. J. Chem. Phys. 2000, 113, 4852–4862. (71) Gelb, L. D. J. Chem. Phys. 2003, 118, 7747–7750. (72) Siepmann, J. I.; Frenkel, D. Mol. Phys. 1992, 75, 59–70. (73) Smit, B.; Karaborni, S.; Siepmann, J. I. J. Chem. Phys. 1995, 102, 2126–2140. (74) Monte Carlo for Complex Chemical Systems-MN. http://chem-siepmann.oit.umn. edu/siepmann/software.html, Accessed: 01-01-2014. (75) Smit, B.; Frenkel, D. Mol. Phys. 1989, 68, 931–950. (76) Rowlinson, J. S.; Widom, B. Molecular Theory of Capillarity; Oxford University Press: New York, 1989; p 261. (77) Rowlinson, J. S.; FL, S. Liquids and Liquid Mixtures; Butterworth: London, 1982; pp 59–85. (78) Rai, N.; Siepmann, J. I. J. Phys. Chem. B 2012, 117, 273–288. (79) Rai, N.; Maginn, E. J. J. Phys. Chem. Lett. 2011, 2, 1439–1443. (80) Rai, N.; Maginn, E. J. Faraday Discuss. 2012, 154, 53–69. (81) Coester, F. Nucl. Phys. 1958, 7, 421–424. (82) Coester, F.; K¨ ummel, H. Nucl. Phys. 1960, 17, 477–485. ˇ ıˇzek, J. J. Chem. Phys. 1966, 45, 4256–4266. (83) C´ 29



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

ˇ ıˇzek, J. Adv. Chem. Phys. 1969, 14, 35. (84) C´ ˇ ıˇzek, J.; Paldus, J. Int. J. Quantum Chem. 1971, 5, 359–379. (85) C´ (86) Raghavachari, K.; Trucks, G. W.; Pople, J. A.; Head-Gordon, M. Chem. Phys. Lett. 1989, 157, 479–483. (87) Dunning Jr, T. H. J. Chem. Phys. 1989, 90, 1007–1023. (88) Kendall, R. A.; Dunning Jr, T. H.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796– 6806. (89) Møller, C.; Plesset, M. S. Phys. Rev. 1934, 46, 618–622. (90) Head-Gordon, M.; Pople, J. A.; Frisch, M. J. Chem. Phys. Lett. 1988, 153, 503–506. (91) Boys, S. F.; Bernardi, F. d. Mol. Phys. 1970, 19, 553–566. (92) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, .; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian 09 Revision C.1. Gaussian Inc. Wallingford CT 2009. 30


Page 30 of 33

Page 31 of 33

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


(93) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Sch¨ utz, M.; Celani, P.; Györffy, W.; Kats, D.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; Shamasundar, K. R.; Adler, T. B.; Amos, R. D.; Bernhardsson, A.; Berning, A.; Cooper, D. L.; Deegan, M. J. O.; Dobbyn, A. J.; Eckert, F.; Goll, E.; Hampel, C.; Hesselmann, A.; Hetzer, G.; Hrenar, T.; Jansen, G.; Köppl, C.; Liu, Y.; Lloyd, A. W.; Mata, R. A.; May, A. J.; McNicholas, S. J.; Meyer, W.; Mura, M. E.; Nicklass, A.; O’Neill, D. P.; Palmieri, P.; Peng, D.; Pfl¨ uger, K.; Pitzer, R.; Reiher, M.; Shiozaki, T.; Stoll, H.; Stone, A. J.; Tarroni, R.; Thorsteinsson, T.; Wang, M. MOLPRO, version 2015.1, a Package of Ab Initio Programs. 2015; see http://www.molpro.net. (94) Nelson Jr, R. D.; Lide Jr, D. R.; Maryott, A. A. Selected Values of Electric Dipole Moments for Molecules in the Gas Phase, National Standard Reference Data SeriesNational Bureu of Standards 10 ; 1967; pp 15–34. (95) Tran, F.; Hutter, J. J. Chem. Phys. 2013, 138, 204103. (96) Lowden, L. J.; Chandler, D. J. Chem. Phys. 1974, 61, 5228–5241. (97) NIST Chemistry Webbook (http://webbook.nist.gov). Accessed: 06-01-2015. (98) Pitzer, K. S. J. Chem. Phys. 1939, 7, 583–590. (99) Prausnitz, J. M.; Lichtenthaler, R. N.; de Azevedo, E. G. Molecular thermodynamics of fluid-phase equilibria; Pearson Education, 1998; pp 104–110. (100) Stone, A. The Theory of Intermolecular Forces; Oxford University Press, 2013; pp 141–148. (101) Björkman, T.; Gulans, A.; Krasheninnikov, A. V.; Nieminen, R. M. Phys. Rev. Lett. 2012, 108, 235502. (102) Björkman, T. Phys. Rev. B 2012, 86, 165109.

31



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

(103) Palmer, B. J.; Anchell, J. L. J. Phys. Chem. 1995, 99, 12239–12248. (104) Hloucha, M.; Deiters, U. K. Fluid Phase Equilibr. 1998, 149, 41–56. (105) Higashi, B. S.; Takada, A. Mol. Phys. 1997, 92, 641–650. (106) Jedlovszky, P.; Mezei, M. J. Chem. Phys. 1999, 110, 2991–3002. (107) Potter, S.; Tildesley, D.; Burgess, A.; Rogers, S. Mol. Phys. 1997, 92, 825–834. (108) Watkins, E. K.; Jorgensen, W. L. J. Phys. Chem. A 2001, 105, 4118–4125. (109) Rai, N.; Siepmann, J. I. J. Phys. Chem. B 2007, 111, 10790–10799.

32


Page 32 of 33

Page 33 of 33

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


Graphical TOC Entry

[tb]

33


Vapor Liquid Equilibria of Hydrofluorocarbons Using Dispersion

Recommend Documents