Optimizing Noble Gas–Water Interactions via Monte Carlo Simulations

Oct 9, 2015 - Each phase consists of 512 molecules (i.e., Ntotal = 1024 molecules) for temperatures ranging between 293.15 and 353.15 K at 1 and 10 at...
0 downloads 10 Views 1MB Size
Subscriber access provided by EPFL | Scientific Information and Libraries

Article

Optimizing Noble Gas-Water Interactions Via Monte Carlo Simulations Oliver Warr, Christopher J. Ballentine, Junju Mu, and Andrew Masters J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.5b06389 • Publication Date (Web): 09 Oct 2015 Downloaded from http://pubs.acs.org on October 10, 2015

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.

The Journal of Physical Chemistry B 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

The Journal of Physical Chemistry

Optimizing Noble Gas-Water Interactions Via Monte Carlo Simulations Oliver Warr1,2*, Chris J. Ballentine1,2, Junju Mu3 and Andrew Masters3

1

School of Earth, Atmospheric and Environmental Sciences, Williamson Building, University

of Manchester, Manchester, M13 9PL, United Kingdom. 2

Department of Earth Sciences, University of Oxford, South Parks Road, Oxford OX1 3AN,

United Kingdom. 3

School of Chemical Engineering and Analytical Science, University of Manchester, Oxford

Road, Manchester M13 9PL, United Kingdom. *

Email: [email protected], phone: +44 1865 272045

KEYWORDS Lennard-Jones optimization, partitioning, GEMC, noble gases

ABSTRACT

1 ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

In this work we present optimized noble gas-water Lennard-Jones 6-12 pair potentials for each noble gas. Given the significantly different atomic nature of water and the noble gases the standard Lorentz-Berthelot mixing rules produce inaccurate unlike molecular interactions between these two species. Consequently we find simulated Henry’s coefficients deviate significantly from their experimental counterparts for the investigated thermodynamic range (293-353 K at 1 and 10 atmospheres), due to a poor unlike potential well term (εij). Where εij is too high or low, so too is the strength of the resultant noble gas-water interaction. This observed inadequacy in using the Lorentz-Berthelot mixing rules is countered in this work by scaling εij for helium, neon, argon and krypton by factors of 0.91, 0.8, 1.1 & 1.05 respectively to reach a much improved agreement with experimental Henry’s coefficients. Due to the highly sensitive nature of the xenon εij term, coupled with the reasonable agreement of the initial values, no scaling factor is applied for this noble gas. These resulting optimized pair potentials also accurately predict partitioning within a CO2-H2O binary phase system as well as diffusion coefficients in ambient water. This further supports the quality of these interaction potentials. Consequently they can now form a well-grounded basis for the future molecular modelling of multiphase geological systems.

1. INTRODUCTION

When trace amounts of noble gases are present within a gas-water binary phase system, they exist as solute particles within both phases1-3. The exact ratio of noble gases in each phase, the partitioning, is dependent on environmental factors such as temperature and phase composition2-3. This relationship can be quantified using Henry’s law: (1) ܲ௜ = ‫ܭ‬௜ሺ்ሻ ‫ݔ‬௜

2 ACS Paragon Plus Environment

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

The Journal of Physical Chemistry

where Pi is the partial pressure of noble gas i in the gaseous phase, Ki(T) is Henry’s coefficient at temperature T and xi is the molar fraction of noble gas in the water phase. This partitioning is used to quantify and constrain fluid interactions for geological systems. They provide information on the formation and evolution of natural geological systems4-5, constrain groundwater recharge rates6-7, as well as allowing us to reconstruct palaeotemperatures

8-9

and model ocean circulation patterns10-11. The partitioning of noble gases can also play an important economic role as it can quantify the origin and evolution of hydrocarbon reservoirs 12-13

and valuable ore deposits14-15.

Their utility in this role is intrinsically linked to an intimate understanding of noble gas behavior within multiphase fluid systems. This can be experimentally calibrated for simple systems but the more extreme complex multiphase systems (e.g. gas-oil-water2 or supercritical CO2-water1) require additional quantitative assessment. Although this can be achieved experimentally1, 16 there are significant drawbacks associated with solely relying on experimentally pre-determined P-V-T-X ranges. For example, one major limitation of this approach is the uncertainty associated with extrapolation between and beyond defined experimental conditions for which noble gas data has been generated (e.g. 2, 17-18).

Due to advances in computing the possibility now exists for efficiently modelling noble gas behavior for specific thermodynamic conditions of interest via molecular simulation19-22. Such an approach therefore sidesteps the typical limitations associated with the experimental approach as a model can be adapted to simulate the exact thermodynamic conditions of interest. However, a robust model requires all intermolecular interactions to be well constrained. Typically for molecular modelling, Lennard-Jones 6-12 interactions between unlike species are obtained using the Lorentz-Berthelot mixing rules23. However, these

3 ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

simplistic calculations have a proven degree of inaccuracy24-27. Here we assess their validity for deriving noble gas-water interactions by using them to simulate Henry’s coefficients for each noble gas over a range of temperatures and pressures (293.15-353.15 K at 1 and 10 atmospheres). Thermally this fits within the original experimental range used to define experimental Henry’s coefficients while the pressures result in fugacity coefficients of 5% or less for the noble gas phase2,

18, 28-29

. Consequently all non-ideal terms (i.e. fugacity and

activity coefficients) are minor and allow a direct comparison with the experimentally derived values. Where discrepancy is observed from experimental values we apply a scaling factor to the εii term of the unlike Lennard-Jones 6-12 interaction to amplify/reduce the noble gas-water interaction. As a result, noble gas-water interactions are iteratively improved using this new, yet simple, technique to generate a series of optimized noble gas-water pair potentials which better replicate experimental observation. As further support for the new proposed terms, additional simulations were run comparing the optimized and original values for a simple carbon dioxide-water binary phase system using a modified Widom insertion technique to derive partition coefficients from the residual chemical potential. Lastly simulated diffusion coefficients, obtained using the new pair potentials, show better agreement with available experimental data than those calculated using the original pair potentials. Consequently these new pair potentials can be used for future simulation with confidence.

2. THEORETICAL METHODS

2.1 Overview Using the isothermal-isobaric (NPT) Gibbs Ensemble Monte Carlo (GEMC) technique30 the Towhee simulation package31 is adapted to simulate a 2-box binary phase system comprised

4 ACS Paragon Plus Environment

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

The Journal of Physical Chemistry

of a water-rich phase (box 1) with an adjacent noble gas phase (box 2). From a starting configuration the system is perturbed in one of three ways: volumes are exchanged between boxes, molecules are moved within each phase and molecules are moved between the two phases30. Both the move and magnitude are randomly determined via a Markov chain random number generator. After the perturbation is implemented the overall change to the distribution of potential energy within the system is evaluated. The new configuration is either accepted or rejected based on the Metropolis algorithm32 to ensure a correct Boltzmann distribution is maintained for the system. To ensure a sufficiently dynamic system a feedback algorithm is used within Towhee to maintain an overall target acceptance ratio for each move type, i.e. the percentage of successful perturbations, which is specified by the user. This is implemented by initially defining the maximum displacement for all applicable perturbations. During a run the cumulative acceptance of each move type (Pacc) is calculated. Where the cumulative acceptance is lower than that specified (Pacc < Target Acceptance) the associated maximum displacement is reduced to increase the probability of acceptance while high acceptances (Pacc > Target Acceptance) result in increased maximum displacement. Additional references are provided to the reader for a full, in-depth description of Towhee, NPT simulations and the GEMC approach23, 30-31, 33. The probability and target acceptance of each type of GEMC move being applied to the system in this study is outlined in Table 1.

Table 1. List of GEMC system moves System move

Probability

Target acceptance

Volumetric

1%a, 1%b

50%a,b

Interphase

62%a, 39%b

N/A

Intraphase

37% a, 60%b

50% a,b

5 ACS Paragon Plus Environment

The Journal of Physical Chemistry

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

List of GEMC system moves with probability of move being accepted and target acceptance where appropriate. ‘Volumetric’ refers to perturbations involving exchange of volume between boxes, ‘interphase’ refers to movement of molecules between boxes and ‘intraphase’ refers to movement of a randomly selected molecule within a box chosen at random. Target acceptance denotes the optimal percentage of moves accepted where applicable and is controlled by a subroutine within Towhee. Values designated ‘a’ and ‘b’ refer to the pure noble gas-water and carbon dioxide-water systems respectively.

Each phase consists of 512 molecules (i.e. Ntotal = 1024 molecules) for temperatures ranging between 293.15-353.15 K at 1 and 10 atmospheres (1.01325 and 10.1325 bar) respectively for each of the five noble gases. All Towhee simulations are run on single processor jobs on the Computational Shared Facility (CSF) cluster at the University of Manchester. To generate a viable starting configuration for deriving Henry’s coefficients each simulation is preequilibrated to allow phase equilibrium for the given thermodynamic conditions. This is achieved by allowing each system to equilibrate for a minimum of 60 million moves. Phase equilibration is assessed by monitoring phase density; when this ceases to fluctuate significantly (