Water System in Predicting

Article Cite This: J. Chem. Eng. Data 2019, 64, 2464−2474

pubs.acs.org/jced

Model Comparison of the CH4/CO2/Water System in Predicting Dynamic and Interfacial Properties Cong Chen,†,‡ Wenfeng Hu,†,‡ Weizhong Li,†,‡ and Yongchen Song*,†,‡ †

Key Laboratory of Ocean Energy Utilization and Energy Conservation of Ministry of Education and ‡School of Energy and Power Engineering, Dalian University of Technology, Dalian 116024, P. R. China

Downloaded via BUFFALO STATE on August 2, 2019 at 22:45:06 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

S Supporting Information *

ABSTRACT: The CO2 storage in deep saline aquifers and enhanced oil/gas recovery techniques have an increasing demand of the knowledge about dynamic properties of natural gas components such as methane and carbon dioxide, as well as the interfacial properties of multiphase fluids such as CH4/water and CO2/water systems. Molecular dynamics simulation is a very important and efficient method to predict dynamic and interfacial properties. The force-field parameters are most important for the accuracy of molecular dynamics simulation. Although there are many potential energy models for CO2/CH4/water, the single model often has high accuracy only for a specific characteristic. It is necessary to find or develop a combined potential energy model that can describe the set of properties with high accuracy simultaneously. In this study, a comprehensive comparison of the quality of different CO2/CH4/water force fields in predicting dynamic and interfacial properties was conducted. The models we selected here included seven different CO2 models (three rigid models and four flexible models), two CH4 models (a united atom (UA) model and an all-atom (AA) model), and three flexible water models. It was found that the rigid CO2 model performs better than its flexible version in density but has little difference in self-diffusion. The UA model of CH4 slightly outperformed the AA model both in density and self-diffusion, and the UA model can save computing resources. The accuracy of interfacial properties depends more on which water model was adopted, and it shows that the F3C water model is the best choice within the water models selected. These results might be helpful to select suitable force fields to further investigate the dynamic and interfacial properties of CO2/CH4/water systems.

1. INTRODUCTION

CO2 is usually injected into geological formation at the depth of more than 800 m to achieve a supercritical state (T > 304.1 K, P > 7.376 MPa). The main mechanisms of CO2 capture in geological structures include structure capture, residual capture, dissolution capture, and mineralization capture.10 The supercritical CO2 will move upward gradually under the force of the buoyancy, and the cap rock which is the first barrier will prevent the upward movement, achieving structure capture, avoiding the leakage of CO2 into the atmosphere. During the migration of CO2 in the reservoir, part of the CO2 is permanently trapped in the pores of the rock because of capillary pressure, achieving residual capture.11 The

In recent years, the excessive use of fossil fuels has caused an energy crisis and brought a series of environmental problem such as the greenhouse effect which seriously threatens the ecological environment.1 The carbon capture and geological sequestration technology has attracted increasing attention, becoming a major strategy to cope with climate change.2 Geological storage sites mainly include: deep saline aquifers, depleted oil and gas reservoirs, infertile coal layer, and geothermal reservoirs.3 The deep saline aquifers are considered as the most potential site because of the wide distribution and the large storage capacity.4,5 In addition, the collected CO2 can be injected into the oil/gas reservoirs to enhance oil/gas recovery (EOR/EGR), achieving permanent geologic storage of CO2 while improving oil/gas recovery.6−9 © 2019 American Chemical Society

Received: January 3, 2019 Accepted: May 23, 2019 Published: June 4, 2019 2464

DOI: 10.1021/acs.jced.9b00006 J. Chem. Eng. Data 2019, 64, 2464−2474

Journal of Chemical & Engineering Data

Article

The interfacial properties of CO2/water31−35 and CH4/ water36−39 systems have been explored using MD simulations with different combined models. Molecular simulation has been regarded as an effective tool to study CO2/water and CH4/water interactions but lack of well-tested CO2/H2O and CH4/H2O force fields.40−46 Therefore, it is necessary to find or develop a combined potential energy model that can describe the characteristics of dynamic and interfacial properties of the CO2/CH4/water system with high accuracy simultaneously. We report here a comparison of CO2 and CH4 force fields, in order to evaluate which force field for CO2 and CH4 performs best in the calculation of density and self-diffusion. At the same time, we have systematically evaluated the force fields of water to simulate the interfacial properties of CO2/water and CH4/ water systems. We aim to provide a certain evidence to choose the molecular force field for simultaneously studying the characteristics of density, self-diffusion coefficient, and IFT.

migration upward includes two mechanisms: molecular diffusion and Darcy flow.12,13 Molecular diffusion occurs everywhere, and Darcy flow occurs when the pressure difference between the wetting phase (water) and no-wetting phase (CO2) is above the capillary pressure of the pore throat of arbitrarily large sizes, in that case, CO2 likely leaks to the atmosphere. Therefore, the capillary pressure is an important parameter to affect the storage security and capacity, and it can be calculated by the standard Young−Laplace14,15 equation Pc = Pg − Pw =

2γwg cos Ø R

(1)

where Pc is the capillary pressure, Pg is the pressure of gas phase, Pw is the pressure of water, γwg is the interfacial tension (IFT) between water and gas, Ø is the wetting angle, and R is the radius of connected pore throat. It shows that the IFT is an important parameter to evaluate whether storage is safe. In the injection of CO2 to EOR/EGR, the IFT also plays a vital role in the storage security and capacity. When CO2 is injected to replace the oil or gas, the IFT of the CO2/water system is much lower than the IFT of the CH4/water system under the same conditions, causing the reduction of capillary pressure of the cap rock. The injected CO2 may migrate into the upper formation by Darcy flow, causing the reduction of storage capacity.13 In order to better carry out the CO2 capture and storage and enhanced oil/gas projects, it is necessary to increase the understanding of the dynamic and interface properties of multiphase fluids in closed structures. In the case of the absence of experimental data in extreme experimental conditions, we can adopt molecular dynamics (MD) simulation technique, a powerful and widely used tool, to predict the physical and chemical properties of fluids but limited by the accuracy of molecular force fields available. There is a large number of force fields proposed for CO2 and CH4,16−22 while the single model often has high accuracy only for a specific characteristic. When investigating the dynamic and interfacial properties of CO2/water and CH4/water systems, we need to simultaneously study the characteristics of density, self-diffusion coefficient, and IFT. Many researchers18,23−26 used MD methods to predict the self-diffusion coefficient and density of pure CO2 with different CO2 models. In addition, the comparison of different models has been conducted by several studies. Zhang27 made a comparison of EPM2 and Zhang models in the prediction of self-diffusion for temperature between 243.0 and 373.0 K and pressures up to 50 MPa, it showed that the Zhang model gives the best prediction of self-diffusion at high temperatures while EPM2 shows better accuracies at 243.0 K. Stassen28 made a comparison of different CH4 potential models for liquid methane at a temperature of 150 K and a density of 28 mol/L, and the result showed that the united atom (UA) model may provide accurate predictions in dynamic and transport properties. Aimoli29 used seven different CO2 models and two different CH4 models to calculate the density of CO2 and CH4 under supercritical conditions up to 900 K and 100 MPa using MD simulations. In a later study of Aimoli,30 they made a comparison of seven different CO2 models and three different CH4 models in predicting the self-diffusion for temperatures ranging from 273.15 to 573.15 K and pressures up to 800 MPa. Different models have good accuracy only at specific temperature and pressure ranges.

2. METHODOLOGY In this work, we first employ different CH4 and CO2 molecular models to calculate the density and self-diffusion of pure fluids at same temperatures and different pressures. Then, we combined CH4 and CO2 models with flexible H2O models, respectively, to calculate the IFT of CH4/water and CO2/water systems. The simulated results are compared to experimental data. 2.1. Molecular Model Used. In the paper, seven CO2 molecular models have been considered, including three threesite rigid models (EPM2,17 TraPPE,47 Zhang27) and four three-site flexible models (Zhu,18 Cygan,48 EPM2-flex,25 TraPPE-flex49). All the models we choose have a better accuracy for the description of a certain characteristic of CO2. The widely used model EPM217 is a very successful rigid model in predicting the liquid−vapor coexistence. EPM217 was modified to a fully flexible version (EPM2-flex25) that included the vibrational modes of the CO bond distance by NietoDraghi, and the EPM2-flex25 model can be used to study the effect of the flexibility of molecular bonds on the thermodynamic properties of CO2. The TraPPE47 model, which is a nonpolarizable transferable potential for phase equilibrium force field, can be used to predict the phase behavior of complex multicomponent systems with a high accuracy. It was modified to a fully flexible model TraPPEflex49 by Perez-Blanco. The Zhang27 model showed advantages in volumetric properties, phase behaviors, and dynamic properties comparing with MSM20,50 and the model proposed by Errington.19 The Zhu18 model has some advantages in reproducing the supercritical property. Cygan48 expanded the previous three-site potential and developed a flexible model to examine the effect of molecular flexibility on the diffusion rate of CO2 in water. There are two different methods to describe the molecular force field of CH4. The first way is to treat each hydrogen and carbon atom as an interaction site51 (AA model), and the second way is to unite each carbon and its bonded hydrogen into a single-interaction site (UA model). Stassen28 presented a comparison of different force fields for liquid methane, and they believed that AA models can provide a more realistic structure properties of liquid methane. However, the UA52 model can give accurate prediction of thermodynamic properties with less computational burden. The vapor−liquid coexistence curve and the critical properties of CH4 calculated by the TraPPE22 model based on the UA description are in 2465



Article

Table 1. Force-Field Parameters for CO2 Models Used models

EPM2-flex

Zhu

Cygan

TraPPE-flex

EPM2

Zhang

TraPPE

ϵC/kcal·mol−1 εO/kcal·mol−1 σC/Å σO/Å lC−O/Å θ0/(deg) qC/e qO/e Kr/kcal mol−1 Å−2 Kθ/kcal mol−1 rad−2

0.055855 0.159860 2.757 3.033 1.149 180 0.6512 0.3256 1282.462 147.5997

0.055876 0.159558 2.800 3.028 1.162 180 0.6516 0.3258 71.65 13.136

0.055884 0.159604 2.800 3.028 1.162 180 0.6512 0.3256 1008.24 53.965

0.053613 0.156868 2.800 3.05 1.160 180 0.7000 0.3500 1028.266 55.96

0.055855 0.159860 2.757 3.033 1.149 180 0.6512 0.3256

0.057277 0.164128 2.7918 3.000 1.163 180 0.5888 0.2944

0.053613 0.156868 2.800 3.05 1.160 180 0.7000 0.3500

excellent agreement with experiment. According to the study of Aimoli,29 the simulation results of UA models (TraPPE22 and SAFT-γ21 force fields) are more accurate than the prediction of PR EoS, and the SAFT-γ model gives underestimate densities as the pressure reaches to 50 MPa, they attributed this difference to SAFT-γ having a slightly larger repulsive exponent which can lead to repulsive forces between the molecules. In this work, we selected the UA model−TraPPE22 model and the fully flexible AA model− OPLS53 model to study the effect of number of sites on the description of dynamic and interfacial properties of CH4/water system. In order to calculate the IFT of CH4/water and CO2/water system, we selected three three-site force fields of water: SPCflex,54 F3C,55 and TIP3P.56 The SPC-flex54 model is a flexible version of the SPC model originally proposed by Teleman.50 The flexible three-centered water model (F3C55) is well-suited to high-speed computation of long MD trajectories of macromolecules in solution because it works with shortrange truncation. The TIP3P56 model is reparameterized the three-site model to improve calculation accuracy of the energy and density for liquid water by Jorgensen. The selection of flexible water models is based on three considerations: (i) the proper treatment of internal vibrations of the flexible model is essential to describe transport properties of polar liquids;54 (ii) the bond flexibility plays a key role in the simulation of water surface tension;57,58 and (iii) the IFT values of the CO2/water system (obtained by rigid models) from Nielsen59 deviates significantly from the experimental value, indicating that rigid models are not suitable for the prediction of IFT. Although there is also the well-known TIP4P/2005 water model proposed by Abascal and Vega,60 this paper only considers the three-point model of water. 2.2. Force Field. We used a CHARMM style force field with the form of the follow eqs 2−7 Utotal = Ucoul + Uvdw + Ustretch + Ubend

Uvdw

The intramolecular interactions contain bond stretching and bond angle bending which are described by harmonic potential.61 The bond stretching potential is given by eq 6 Ustretch = K r(r − r0)2

(6)

and the bond angle bending potential is given by eq 7 Ubend = K θ(θ − θ0)2

(7)

where r, r0, θ, and θ0 represent the measured bond length, the equilibrium bond length, the measured angle, and the equilibrium angle, respectively, and Kr and Kθ are the bond stretching and angle bending force constants, respectively. The parameters of different force fields are listed in Tables 1−3. Table 2. Force-Field Parameters for H2O Models Used models −1

ϵH/kcal·mol εO/kcal·mol−1 σH/Å σO/Å lH−O/Å θ0/(deg) qH/e qO/e Kr/kcal mol−1 Å−2 Kθ/kcal mol−1 rad−2

(2)

SPC-flex

F3C

TIP3P

0 0.1554 0 3.166 0.96 108.47 0.4238 0.8476 540.6336 50

0.01 0.1848 0.820 3.166 1.0 109.47 0.41 0.82 250 60

0.046 0.1521 0.4000 3.1507 0.9572 104.52 0.417 0.834 450 55

2.3. Simulation Details. All simulations were performed with the LAMMPS62 package. The periodic boundary conditions were applied in three dimensions in all simulations. The initial velocity distribution of molecules satisfied the Maxwell−Boltzmann distribution. The nonbonded interaction was described by L-J potential function. The cut-off radius was set to be 1.4 nm and less than half of the minimum of the three-dimensional size. A particle−particle−particle−mesh63 method and the desired relative error in forces of 0.0001 were used to calculate the long-range Coulombic interactions. The

qiqj 4πε0rij

(4)

where εij is the well depth for short-range interactions and σij is the core diameter for L-J potential. The L-J parameters for different atoms are determined by Lorentz-Berthelot combining rules σi + σj σij = εij = εiεj (5) 2

The Coulombic interactions are calculated according to Coulomb’s law Ucoul =

ÅÄÅ ÑÉ ÅÅij σ yz12 ij σ yz6ÑÑÑ ij z Ñ ÅÅjj ij zz j = 4εijÅÅÅjj zz − jjj zzz ÑÑÑÑ j rij z ÑÑ ÅÅj rij z k { ÑÑÖ ÅÅÇk {

(3)

where qi and qj are the partial charges of atom i and j, ε0 is the permittivity of free space, and rij is the distance between atoms. The van der Waals interactions, which are described by Lennard-Jones (L-J) potential 2466



Article

with CO2 in the middle and water on both sides was built, where the box length in the z direction was made 50% larger than that in x and y directions. The number of water molecules was fixed to 2304, and the dimension of the water box is 40 Å × 40 Å × 43.5 Å; the number of CO2 molecules was always adjusted according to the pressures. Therefore, the dimension of the whole box was 40 Å × 40 Å × Z, the dimension of the Z direction was determined according to the pressures. A total of 11 ns run was performed, including the first 9 ns in the NPT ensemble to equilibrate and the last 2 ns in the NVT ensemble to produce data. The calculation process for CH4 is similar to that of CO2. We compared the models by calculating the absolute deviations (AD) or average AD (AAD) between the simulated and experimental values via the eqs 8 and 9

Table 3. Force-Field Parameters for CH4 Models Used UA (TraPPE)

AA (OPLS)

ϵCH4/kcal·mol−1

models

0.294076

0.066

εH/kcal·mol−1 σCH4/Å

3.73

0.03 3.5

σH/Å lC−H/Å θ0/(deg) qC/e qH/e Kr/kcal mol−1 Å−2 Kθ/kcal mol−1 rad−2

2.5 1.09 107.8 −0.24 0.06 340 33

temperature and pressure were fixed using Nosé−Hoover64 thermostat and barostat. The time step was 1 fs, the time step to calculate the nonbonded interactions was 1 fs, and the time step to calculate the Coulombic interactions was 2 fs. The molecular graphics software Visual Molecular Dynamics65 was used to visualize and analyze simulation results using in-house codes. For the calculation of CO2 density and self-diffusion coefficient, the number of CO2 molecule was fixed to 3136, the initial dimension of the cubic box simulation was 150 Å × 150 Å × 150 Å. First, energy minimization had been done prior to the NPT simulations, then we ran 2 ns in the NPT ensemble to make the system reach equilibrium. The configuration after equilibrium corresponded to a minimum energy. The dimension of the cubic box after equilibrium was different under different pressures, and it was followed by a production run of 0.5 ns, during which the trajectories were stored every 0.5 ps. Therefore, a total of 1000 trajectory files were generated and then were split into 5 equal intervals, and the standard deviations of the mean were used to estimate the uncertainties. In order to calculate the self-diffusion coefficients of CO2, the equilibrium configuration from NPT simulation was used as an input to perform MD simulations in the NVE ensemble for 1 ns. The trajectories were stored every 1 ps, and then a total of 1000 trajectory files were split into 5 blocks to predict mean and standard deviations of the mean. For the calculation of IFT of the CO2/water system, a rectangular box

AD =

Aiexp − Aical 100% Aiexp

(8)

N

AAD =

exp cal 1 i Ai − Ai 100% ∑ Ni i = 1 Aiexp

(9)

where AD is the absolute deviation, AAD is the average absolute deviation, Aexp is the experimental values, and Acal i i is the simulated values.

3. RESULTS AND DISCUSSION 3.1. Density. Density is one of the most important physical properties of fluids and here was calculated via the eq 10 ρ=

NM ⟨V ⟩

(10)

where N is the number of CO2 or CH4 molecules in the system, M is the relative molecule mass of CO2 or CH4, V is the volume of the system, and ⟨...⟩ means ensemble average. MD simulations in the NPT ensemble were adopted to calculate the density of pure CO2 and CH4 at a temperature of 323 K and pressure ranged between 2 and 21 MPa. We compared the simulated results of all models with NIST data.

Figure 1. Density isotherms corresponding to T = 323 K of different CO2 models (left) and the AD between simulated and experimental data (right). The error bars for simulated densities are smaller than the symbol size. 2467



Article

Figure 2. Density isotherms corresponding to T = 323 K of different CH4 models (left) and the AD between simulated and experimental data (right). The error bars for simulated densities are smaller than the symbol size.

Table 4. Overall AAD of CO2 and CH4 Properties Calculated through Molecular Simulations When Compared to Experimental Data AAD (%) density TraPPE Zhang EPM2 Cygan-flex EPM2-flex TraPPE-flex Zhu-flex UA AA

Figure 3. Temperature of the system varies with time during NVE simulations with the EPM2 model at 3.4 MPa.

CO2 7.92 11.74 11.18 13.15 16.01 12.56 10.38 CH4 2.72 5.14

self-diffusion 12.33 19.93 22.13 18.21 22.96 15.31 16.13

The simulated density data of CO2 and CH4 are summarized in Tables S1 and S2 in the Supporting Information.

Figure 4. Self-diffusion isotherms corresponding to T = 323 K of different CO2 models (left) and the AD between simulated and experimental data (right). The error bars for simulated self-diffusion coefficients are smaller than the symbol size. The estimated errors of the experimental results are ±0.9%.66 2468



Article

processing of intramolecular bonds. Figure 1 shows that the TraPPE model is obviously better than its flexible version TraPPE-flex, and EPM2 is also better than EPM2-flex in predicting density, indicating that the bond stretching and the bong angle bending of the molecule will increase the instability of the system, leading to an increase in simulation bias, which is consistent with the research results of Aimoli.29 The flexible Zhu model can give relative accurate simulation results because of its smaller bond constants than other flexible models. The results of CH4 density using UA (TraPPE) and AA (OPLS) models are given in Figure 2. It can be seen that as the pressure increases, the density almost increases linearly. Both force fields underestimate densities, but the results of the UA model are closer to the experimental data than those of the AA model. This may be due to the consideration of short-range force increases the instability of the gas system. In addition, because the UA model does not consider the bond interactions of the molecule, its computational efficiency is about three times that of the AA model. Therefore, the UA model is more suitable to predict the property of density. It should be noted that the above conclusions for CO2 and CH4 density comparison between simulation and experiments apply even considering the uncertainties. 3.2. Self-Diffusion. The self-diffusion coefficient is a vital parameter to describe the fluid transport properties. Here, we used an equivalent approach based the Einstein relation to calculate the self-diffusion. Eq 11 constantly calculates the mean square displacement of CO2 or CH4 molecules over time and then derives the time.

Figure 5. Self-diffusion isotherms corresponding to T = 323 K of different CH4 models, and the experimental data of experiment 1 are from Harris67 and of experiment 2 are from Dawson.68 The error bars for simulated self-diffusion coefficients are smaller than the symbol size. The estimated errors of the experimental results are ±267 and ±6%.68

Figure 1 gives the CO2 density results obtained from three rigid models (TraPPE, Zhang, EPM2) and four flexible models (Cygan-flex, EPM2-flex, TraPPE-flex, Zhu). All models can give accuracy prediction of density below the critical pressure of 7.376 MPa, but when the pressure increases above the critical pressure, the simulated values of the density are lower than the experimental values. The deviation between simulated and experimental values increases gradually and it presents the largest deviation at 10.14 MPa. The TraPPE model yields more accurate results than other models as the pressure increases. The main difference between TraPPE and TraPPEflex models and EPM2 and EPM2-flex models lies in the

D = lim

t →∞

1 6Nmt

Nm

∑ [rj(t ) − rj(0)]2 j=1

(11)

where D is the self-diffusion coefficient, Nm is the number of particles, rj(t) is the position of j particle at time t, t is the simulation time, and ⟨...⟩ means ensemble average. The

Figure 6. Final Snapshot of the simulation box of the CO2/water system using the EPM2-flex−SPC-flex model at 333.15 K and 15 MPa and the density profile along the z axis and the energy change of the system (left). The energy change of the system (right). 2469



Article

Figure 7. Pressure dependence of IFT at 333.15 K calculated by MD simulation of (a) Cygan−SPC (flex), TraPPE−SPC (flex), EPM2−SPC (flex), Zhu−SPC (flex), and (b) Cygan−F3C, TraPPE−F3C, EPM2−F3C, Zhu−F3C, and (c) Cygan−TIP3P, TraPPE−TIP3P, EPM2−TIP3P, Zhu−TIP3P for the CO2/water system. The experimental data of experiment 1 are from Hebach and69 experiment 2 are from Georgiadis.70 The error bars for simulated IFTs are smaller than the symbol size. The AAD between simulated and experimental 1 data are listed at bottom.

temperature of the system was fixed at 323 K for all simulations and the pressure ranged from 2 to 21 MPa. Figure 3 presented the temperature of the system varies with time during NVE simulations with the EPM2 model at 3.4 MPa, and the average temperature is 323.43 ± 1.03 K. This shows that the system reached the equilibrium. Other systems used the same method to determine whether the system was equilibrated. We compared the results of all models with experimental values from Etesse66 for CO2 and from Harris67 and Dawson68 for CH4. The simulated self-diffusion coefficient data of CO2 and CH4 are summarized in Tables S3 and S4 in the Supporting Information. The self-diffusion isotherms of CO2 in Figure 4 show that the self-diffusion decreases with the increase of pressure, and the slope of the curve decreases gradually with increasing

pressure. When the pressure exceeds the critical value, the selfdiffusion decreases slowly. All models overestimate the selfdiffusion coefficient; in addition, there is no significant difference in the results between the rigid and flexible CO2 model, indicating that the self-diffusion coefficient is insensitive to the bond and angle flexibility. According to Figure 4 (right), we can see that it also presents the largest deviation at 10.14 MPa. The overall AAD of CO2 density and self-diffusion coefficient calculated through molecular simulations when compared to experimental data are presented in Table 4. It can be seen that the TraPPE model is the most accurate, followed its flexible version TraPPE-flex model, and the flexible Zhu model provides similar results as the TraPPE-flex model. The self-diffusion isotherms (Figure 5) of CH4 have the same trend with the curve of CO2, and the slope of the curve 2470



Article

Figure 8. Pressure dependence of IFT at 333.15 K calculated by the MD simulation of (a) UA−SPC (flex), AA−SPC (flex), UA−TIP3P, AA− TIP3P, and UA−F3C, AA−F3C for the CH4/water system. The experimental data of experiment are from Ren.71 (b) Renormalizing the predicted IFT values provided by CH4−SPC, CH4−TIP3P, and CH4−F3C models to the surface tension of water shifts simulated values upward by 14.69, 11.98, and 11.65 mN/m, respectively. The error bars for simulated IFTs are smaller than the symbol size. The AAD between simulated and experimental data based on renormalized data are listed at bottom.

where Pzz is pressure in the direction normal to the CO2/water or CH4/water interface, Pxx and Pyy are pressures in directions parallel to the interface, and lz is the simulation box length in the direction normal to the interface. The simulated IFT data of the CO2/water system and CH4/water system are summarized in Tables S5 and S6 in the Supporting Information. The simulation results of the CO2/water system in our study are presented in Figure 7a−c, and it shows that the curves of IFT have similar trend with the experimental values. All combined models show that the IFT decreases significantly with pressure below the critical pressure of 7.376 MPa and then decreases slowly with pressure. When we combined different CO2 models with the same water model, the curve of pressure dependence of IFT tends to be consistent, indicating that the curve depends more on which the water model was used instead of the CO2 model. The CO2−SPC (flex) model and CO2−TIP3P model provide high accurate results under high pressures, but the results under low pressures are not satisfactory. The curves yielded by CO 2 −F3C model combinations run approximately parallel to the experiment curve under all pressures, showing that the interfacial property of CO2/water was best captured when the F3C water model is used instead of SPC (flex) and TIP3P water models. According to the AAD, the Zhu−F3C model gives best results among all simulations. The simulation result of the CH4/water system in our study is presented in Figure 8a. The IFT of CH4/water decreases

decreases gradually with the increasing pressure. Because the experimental data and the simulated data are not under the same pressure condition, the AD cannot be obtained. However, it is obvious from Figure 5, it shows that the UA (TraPPE) model provided more accurate prediction of selfdiffusion apparently than the AA (OPLS) model. In the case of a given computing resource, we prefer to choose the UA model to predict the property of diffusion. The overall AAD of CH4 density calculated through molecular simulations when compared to experimental data are presented in Table 4. The overall AAD of CH4 density also prove that the UA model is a better choice than the AA model. 3.3. IFT of CO2/Water and CH4/Water Systems. MD simulations in the NVT ensemble were adopted to calculate the IFT with different force-field models at temperature of 333.15 K and different pressures. The density profiles of water/ CO2 and a snapshot of the equilibrated simulation box are shown in Figure 6. It can be seen that the density fluctuates very little along in the direction normal to the interface, and CO2 has been uniformly distributed in water. In addition, the energy of the system during the data production fluctuates up and down the equilibrium value (Figure 6). These show good equilibration of simulation systems. The value of IFT was calculated by the principle components of MD simulation cell stress tensor Ä ÉÑ l ÅÅ Ñ 1 γ = z ÅÅÅÅPzz − (Pxx + Pyy)ÑÑÑÑ ÑÖ (12) 2 ÅÇ 2 2471



■

almost linearly with the increasing of pressure and all simulated values underestimate the experimental values. The surface tensions of water at 333.15 K predicted by SPC (flex), TIP3P, and F3C water models were calculated and the values were 51.55, 54.59, and 54.26 mN/m, respectively. The surface tension of water reported by experiment was 66.24 mN/m.72 The differences between experimental and simulated values obtained by SPC (flex), TIP3P, and F3C water models are 14.69, 11.98, 11.65 mN/m, respectively. To account the effect of water models on the IFT, we renormalized the predicted IFT of the CH4/water system by adding the differences between experimental and simulated values of surface tensions (Figure 8b). It was found that the deviations between experimental values and simulated values increase with the increasing pressures. AAD showed that the interfacial property of CH4/water was best captured when the AA−F3C model is used, this indicated the importance of proper treating internal vibrations of the flexible model for describing transport properties. However, considering the higher efficiency and the better prediction in density and diffusion of the UA model, we prefer to adopt the UA−F3C model to investigate the properties of CH4/water systems. The simulated IFT data of the CH4/water system are summarized in Table S6 in the Supporting Information.

Article

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jced.9b00006.

■

Simulated densities and self-diffusion coefficients of CH4 and CO2 and the simulated IFTs of CO2/water and CH4/water systems (PDF)

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Cong Chen: 0000-0001-5035-1026 Funding

This research was supported by the National Natural Science Foundation of China (51676027 and 51206016), the Natural Science Foundation of Liaoning Province (201602147), and the Fundamental Research Funds for the Central Universities (DUT14LAB13 and DUT17LAB03). Notes

The authors declare no competing financial interest.

■

4. CONCLUSIONS In the present study, a comprehensive comparison of the quality of different CO2 and CH4 force fields in predicting dynamic properties such as density and self-diffusion coefficient and interfacial properties is conducted. The reliability of different models depends on the conditions studied and the required characteristics. For CO2, the rigid TraPPE model performed best both in density and selfdiffusion. The density given by the rigid models generally outperforms the flexible models, suggesting that intramolecular interactions may increase the instability of the gas system; however, the flexible Zhu model gives relative accurate simulation results of density because of its smaller bond constants. Unlike density, intramolecular interactions have little effect on the self-diffusion coefficient. For CH4, the UA model performs better than the AA model in predicting both density and self-diffusion coefficients, and the computational efficiency of the UA model is higher. For the CO2/water system, the accuracy of IFT depends more on which the water model was used instead of the CO2 model. The interfacial property of CO2/water was best captured when the F3C water model is used, and the Zhu−F3C model gives results closest to the experiment. For the CH4/water system, the IFTs were better predicted when the F3C model was used with considering the deviation of predicted and experimental water surface tension. Therefore, the Zhu−F3C−UA model can give good estimations for density, diffusion, and interfacial characteristics of the CO2/CH4/water system. In CO2 sequestration industry and enhanced oil/gas projects, we need to obtain the dynamic and interfacial properties of the CO2/CH4/water system simultaneously when studying the migration and adsorption characteristics of fluids. We believe that these studies could be helpful to better understand the limits of force fields in describing the dynamic and interfacial properties of CO2/water and CH4/water systems.

REFERENCES

(1) Oreskes, N. The scientific consensus on climate change. Science 2004, 306, 1686. (2) Benson, S. M.; Surles, T. Carbon dioxide capture and storage: An overview with emphasis on capture and storage in deep geological formations. Proc. IEEE 2006, 94, 1795−1805. (3) Bachu, S. From suitability to ultimate capacity: A roadmap for assessing sedimentary basins and selecting sites for CO2 storage in geological media. Greenhouse Gas Control Technologies: Proceeding of Fifth International Conference; CSIRO Publishing: Collingwood, VIC, Australia, 2001; pp 328−333. (4) Bachu, S. Sequestration of CO2 in geological media: criteria and approach for site selection in response to climate change. Energy Convers. Manage. 2000, 41, 953−970. (5) Szulczewski, M. L.; MacMinn, C. W.; Herzog, H. J.; Juanes, R. Lifetime of carbon capture and storage as a climate-change mitigation technology. Proc. Natl. Acad. Sci. U.S.A. 2012, 109, 5185−5189. (6) Leach, A.; Mason, C. F.; Veld, K. v. t. Co-optimization of enhanced oil recovery and carbon sequestration. Resour. Energy Econ. 2011, 33, 893−912. (7) Gaspar Ravagnani, A. T. F. S.; Ligero, E. L.; Suslick, S. B. CO2 sequestration through enhanced oil recovery in a mature oil field. J. Pet. Sci. Eng. 2009, 65, 129−138. (8) Ettehadtavakkol, A.; Lake, L. W.; Bryant, S. L. CO2-EOR and storage design optimization. Int. J. Greenhouse Gas Control 2014, 25, 79−92. (9) van t Veld, K.; Mason, C. F.; Leach, A. The economics of CO2 sequestration through enhanced oil recovery. Energy Procedia 2013, 37, 6909−6919. (10) Iglauer, S. CO2-Water-Rock Wettability: variability, influencing factors, and implications for CO2 geostorage. Acc. Chem. Res. 2017, 50, 1134−1142. (11) Wang, J. G.; Ju, Y.; Gao, F.; Liu, J. A simple approach for the estimation of CO2 penetration depth into a caprock layer. J. Rock Mech. Geotech. Eng. 2016, 8, 75−86. (12) Li, S.; Dong, M.; Li, Z.; Huang, S.; Qing, H.; Nickel, E. Gas breakthrough pressure for hydrocarbon reservoir seal rocks: implications for the security of long-term CO2 storage in the Weyburn field. Geofluids 2005, 5, 326−334. (13) Li, Z.; Dong, M.; Li, S.; Huang, S. CO2 sequestration in depleted oil and gas reservoirs-caprock characterization and storage capacity. Energy Convers. Manage. 2006, 47, 1372−1382. 2472



Article

(35) Li, X.; Ross, D. A.; Trusler, J. P. M.; Maitland, G. C.; Boek, E. S. Molecular Dynamics Simulations of CO2 and Brine Interfacial Tension at High Temperatures and Pressures. J. Phys. Chem. B 2013, 117, 5647−5652. (36) Sakamaki, R.; Sum, A. K.; Narumi, T.; Ohmura, R.; Yasuoka, K. Thermodynamic properties of methane/water interface predicted by molecular dynamics simulations. J. Chem. Phys. 2011, 134, 144702. (37) Reed, S. K.; Westacott, R. E. The interface between water and a hydrophobic gas. Phys. Chem. Chem. Phys. 2008, 10, 4614−4622. (38) Biscay, F.; Ghoufi, A.; Lachet, V.; Malfreyt, P. Monte Carlo Simulations of the Pressure Dependence of the Water−Acid Gas Interfacial Tensions. J. Phys. Chem. B 2009, 113, 14277−14290. (39) Ghoufi, A.; Malfreyt, P. Numerical evidence of the formation of a thin microscopic film of methane at the water surface: a free energy calculation. Phys. Chem. Chem. Phys. 2010, 12, 5203−5205. (40) Vega, C.; de Miguel, E. Surface tension of the most popular models of water by using the test-area simulation method. J. Chem. Phys. 2007, 126, 154707. (41) Lafitte, T.; Mendiboure, B.; Piñeiro, M. M.; Bessières, D.; Miqueu, C. Interfacial Properties of Water/CO2: A Comprehensive Description through a Gradient Theory-SAFT-VR Mie Approach. J. Phys. Chem. B 2010, 114, 11110−11116. (42) Miqueu, C.; Míguez, J. M.; Piñeiro, M. M.; Lafitte, T.; Mendiboure, B. Simultaneous application of the gradient theory and monte carlo molecular simulation for the investigation of methane/ water interfacial properties. J. Phys. Chem. B 2011, 115, 9618−9625. (43) de Lara, L. S.; Michelon, M. F.; Miranda, C. R. Molecular Dynamics Studies of Fluid/Oil Interfaces for Improved Oil Recovery Processes. J. Phys. Chem. B 2012, 116, 14667−14676. (44) Míguez, J. M.; Pineiro, M. M.; Blas, F. J. Influence of the longrange corrections on the interfacial properties of molecular models using Monte Carlo simulation. J. Chem. Phys. 2013, 138, 034707. (45) Míguez, J. M.; Garrido, J. M.; Blas, F. J.; Segura, H.; Mejía, A.; Piñeiro, M. M. Comprehensive characterization of interfacial behavior for the mixture CO2 + H2O + CH4: comparison between atomistic and coarse grained molecular simulation models and density gradient theory. J. Phys. Chem. C 2014, 118, 24504−24519. (46) Müller, E. A.; Mejia, A. Comprehensive resolving discrepancies in the measurements of the interfacial tension for the CO2 + H2O mixture by computer simulation. J. Phys. Chem. Lett. 2014, 5, 1267− 1271. (47) Potoff, J. J.; Siepmann, J. I. Vapor-liquid equilibria of mixtures containing alkanes, carbon dioxide, and nitrogen. AIChE J. 2001, 47, 1676−1682. (48) Myshakin, E. M.; Saidi, W. A.; Romanov, V. N.; Cygan, R. T.; Jordan, K. D. Molecular dynamics simulations of carbon dioxide intercalation in hydrated Na-Montmorillonite. J. Phys. Chem. C 2013, 117, 11028−11039. (49) Perez-Blanco, M. E.; Maginn, E. J. Molecular dynamics simulations of CO2 at an ionic liquid interface: adsorption, ordering, and interfacial crossing. J. Phys. Chem. B 2010, 114, 11827−11837. (50) Murthy, C. S.; Singer, K.; McDonald, I. R. Interaction site models for carbon dioxide. Mol. Phys. 1981, 44, 135−143. (51) Williams, D. E. Nonbonded Potential Parameters Derived from Crystalline Hydrocarbons. J. Chem. Phys. 1967, 47, 4680. (52) Nagy, J.; Weaver, D. F.; Smith, V. H. A Comprehensive Study of Alkane Nonbonded Empirical Force Fields. Suggestions for Improved Parameter Sets. J. Phys. Chem. 1995, 99, 8058−8065. (53) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and testing of the OPLS all-atom force field on conformational energetics and properties of organic liquids. J. Am. Chem. Soc. 1996, 118, 11225−11236. (54) Teleman, O.; Jönsson, B.; Engström, S. A molecular dynamics simulation of a water model with intramolecular degrees of freedom. Mol. Phys. 1987, 60, 193−203. (55) Levitt, M.; Hirshberg, M.; Sharon, R.; Laidig, K. E.; Daggett, V. Calibration and testing of a water model for simulation of the molecular dynamics of proteins and nucleic acids in solution. J. Phys. Chem. B 1997, 101, 5051−5061.

(14) Jung, J.-W.; Wan, J. Supercritical CO2 and Ionic Strength Effects on Wettability of Silica Surfaces: Equilibrium Contact Angle Measurements. Energy Fuels 2012, 26, 6053−6059. (15) Chalbaud, C.; Robin, M.; Lombard, J.-M.; Martin, F.; Egermann, P.; Bertin, H. Interfacial tension measurements and wettability evaluation for geological CO2 storage. Adv. Water Resour. 2009, 32, 98−109. (16) Zhang, Y.; Yang, J.; Yu, Y.-X. Dielectric constant and density dependence of the structure of supercritical carbon dioxide using a new modified empirical potential model: A Monte Carlo simulation study. J. Phys. Chem. B 2005, 109, 13375−13382. (17) Harris, J. G.; Yung, K. H. Carbon Dioxide’s Liquid-Vapor Coexistence Curve And Critical Properties as Predicted by a Simple Molecular Model. J. Phys. Chem. 1995, 99, 12021−12024. (18) Zhu, A.; Zhang, X.; Liu, Q.; Zhang, Q. A fully flexible potential model for carbon dioxide. Chin. J. Chem. Eng. 2009, 17, 268−272. (19) Potoff, J. J.; Errington, J. R.; Panagiotopoulos, A. Z. Molecular simulation of phase equilibria for mixtures of polar and non-polar components. Mol. Phys. 1999, 97, 1073−1083. (20) Murthy, C. S.; O’Shea, S. F.; McDonald, I. R. Electrostatic interactions in molecular crystals. Mol. Phys. 1983, 50, 531−541. (21) Lafitte, T.; Apostolakou, A.; Avendaño, C.; Galindo, A.; Adjiman, C. S.; Müller, E. A.; Jackson, G. Accurate statistical associating fluid theory for chain molecules formed from Mie segments. J. Chem. Phys. 2013, 139, 154504. (22) Martin, M. G.; Siepmann, J. I. Transferable Potentials for Phase Equilibria. 1. United-Atom Description of n-Alkanes. J. Phys. Chem. B 1998, 102, 2569−2577. (23) Higashi, H.; Iwai, Y.; Arai, Y. Calculation of self-diffusion and tracer diffusion coefficients near the critical point of carbon dioxide using molecular dynamics simulation. Ind. Eng. Chem. Res. 2000, 39, 4567−4570. (24) Fernández, G. A.; Vrabec, J.; Hasse, H. Self-diffusion and binary Maxwell-Stefan diffusion coefficients of quadrupolar real fluids from molecular simulation. Int. J. Thermophys. 2005, 26, 1389−1407. (25) Nieto-Draghi, C.; de Bruin, T.; Pérez-Pellitero, J.; Bonet Avalos, J.; Mackie, A. D. Thermodynamic and transport properties of carbon dioxide from molecular simulation. J. Chem. Phys. 2007, 126, 064509. (26) Avendaño, C.; Lafitte, T.; Galindo, A.; Adjiman, C. S.; Jackson, G.; Muller, E. A. SAFT-γ force field for the simulation of molecular fluids. 1. A single-site coarse grained model of carbon dioxide. J. Phys. Chem. B 2011, 115, 11154−11169. (27) Zhang, Z.; Duan, Z. An optimized molecular potential for carbon dioxide. J. Chem. Phys. 2005, 122, 214507. (28) Stassen, H. On the pair potential in dense fluid methane. J. Mol. Struct. 1999, 464, 107−119. (29) Aimoli, C. G.; Maginn, E. J.; Abreu, C. R. A. Force field comparison and thermodynamic property calculation of supercritical CO2 and CH4 using molecular dynamics simulations. Fluid Phase Equilib. 2014, 368, 80−90. (30) Aimoli, C. G.; Maginn, E. J.; Abreu, C. R. A. Transport properties of carbon dioxide and methane from molecular dynamics simulations. J. Chem. Phys. 2014, 141, 134101. (31) da Rocha, S. R. P.; Johnston, K. P.; Westacott, R. E.; Rossky, P. J. Molecular Structure of the Water−Supercritical CO2 Interface. J. Phys. Chem. B 2001, 105, 12092−12104. (32) Kvamme, B.; Kuznetsova, T.; Hebach, A.; Oberhof, A.; Lunde, E. Measurements and modelling of interfacial tension for water + carbon dioxide systems at elevated pressures. Comput. Mater. Sci. 2007, 38, 506−513. (33) Banerjee, S.; Hassenklöver, E.; Kleijn, J. M.; Cohen Stuart, M. A.; Leermakers, F. A. M. Interfacial tension and wettability in watercarbon dioxide systems: experiments and self-consistent field modeling. J. Phys. Chem. B 2013, 117, 8524−8535. (34) Arif, M.; Al-Yaseri, A. Z.; Barifcani, A.; Lebedev, M.; Iglauer, S. Impact of pressure and temperature on CO2-brine-mica contact angles and CO2-brine interfacial tension: Implications for carbon geosequestration. J. Colloid Interface Sci. 2016, 462, 208−215. 2473



Article

(56) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1983, 79, 926−935. (57) Yuet, P. K.; Blankschtein, D. Molecular dynamics simulation study of water surfaces: comparison of flexible water models. J. Phys. Chem. B 2010, 114, 13786−13795. (58) López-Lemus, J.; Chapela, G. A.; Alejandre, J. Effect of flexibility on surface tension and coexisting densities of water. J. Chem. Phys. 2008, 128, 174703. (59) Nielsen, L. C.; Bourg, I. C.; Sposito, G. Predicting CO2-water interfacial tension under pressure and temperature conditions of geologic CO2 storage. Geochim. Cosmochim. Acta 2012, 81, 28−38. (60) Abascal, J. L. F.; Vega, C. A general purpose model for the condensed phases of water: TIP4P/2005. J. Chem. Phys. 2005, 123, 234505. (61) Nath, S. K.; Escobedo, F. A.; de Pablo, J. J.; Patramai, I. Simulation of Vapor−Liquid Equilibria for Alkane Mixtures. Ind. Eng. Chem. Res. 1998, 37, 3195−3202. (62) Plimpton, S. Fast Parallel Algorithms for Short-Range Molecular Dynamics. J. Comput. Phys. 1995, 117, 1−19. (63) Hardy, D. J. Multilevel summation for the fast evaluation of forces for the simulation of biomolecules. Ph.D. Thesis, Univ of Illinois at Urbana−Champaign, 2006. Also Department of Computer Science Report No. UIUCDCS-R-2006-2546, May 2006. (64) Martyna, G. J.; Klein, M. L.; Tuckerman, M. Nosé-Hoover chains: The canonical ensemble via continuous dynamics. J. Chem. Phys. 1992, 97, 2635−2643. (65) Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual molecular dynamics. J. Mol. Graphics 1996, 14, 33−38. (66) Etesse, P.; Zega, J. A.; Kobayashi, R. High-pressure nuclearmagnetic-resonance measurement of spin-lattice relaxation and selfdiffusion in carbon dioxide. J. Chem. Phys. 1992, 97, 2022−2029. (67) Harris, K. R. The density dependence of the self-diffusion coefficient of methane at −50°, 25° and 50 °C. Phys. A 1978, 94, 448−464. (68) Dawson, R.; Khoury, F.; Kobayashi, R. Self-diffusion measurements in methane by pulsed nuclear magnetic resonance. AIChE J. 1970, 16, 725−729. (69) Hebach, A.; Oberhof, A.; Dahmen, N.; Kögel, A.; Ederer, H.; Dinjus, E. Interfacial Tension at Elevated Pressures Measurements and Correlations in the Water + Carbon Dioxide System. J. Chem. Eng. Data 2002, 47, 1540−1546. (70) Georgiadis, A.; Maitland, G.; Trusler, J. P. M.; Bismarck, A. Interfacial tension measurements of the (H2O + CO2) system at elevated pressures and temperatures. J. Chem. Eng. Data 2010, 55, 4168−4175. (71) Ren, Q.-Y.; Chen, G.-J.; Yan, W.; Guo, T.-M. Interfacial tension of (CO2 + CH4) plus water from 298 K to 373 K and pressures up to 30 MPa. J. Chem. Eng. Data 2000, 45, 610−612. (72) Vargaftik, N. B.; Volkov, B. N.; Voljak, L. D. International Tables of the Surface Tension of Water. J. Phys. Chem. Ref. Data 1983, 12, 817−820.

2474


Water System in Predicting

Recommend Documents